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Abstract 

Solar  flares  release  energy,  primarily  at  X-ray  and  EUV  wavelengths,  which  is 
then  absorbed  mainly  in  Earth’s  ionosphere.  This  non-uniform  absorption  of  energy 
alters  the  ionosphere’s  structure  and  can  change  the  propagation  of  electromagnetic 
waves  causing  errors  in  GPS  navigation,  false  radar  echoes,  and  loss  of  High  Frequency 
(HF)  radio  communications.  Accurately  modeling  the  ionospheric  response  to  flares 
is  the  Erst  step  in  predicting,  and  then  mitigating,  their  effects.  Accurately  modeling 
these  effects  requires  solar  irradiance  at  a  high  cadence,  to  capture  the  flare,  which 
may  only  last  minutes,  as  well  as  treatment  of  the  effects  of  photoclectrons.  Sami2  is 
Another  Model  of  the  Ionosphere  (SAMI2)  is  a  two-dimensional  ionospheric  model. 
In  its  standard  form,  it  is  unsuitable  for  studying  solar  flare  effects  because  it  relies 
on  a  daily  proxy  to  specify  the  solar  irradiance.  The  model  was  successfully  modi¬ 
fied  to  include:  secondary  ionization  of  photoelectrons  as  well  as  high-cadence  solar 
irradiance  from  the  Flare  Irradiance  Spectral  Model  (FISM),  and  a  new  photoelec¬ 
tron  heating  parameterization.  The  modified  SAMI2  was  successfully  used  to  model 
the  ionospheric  response  to  the  X17  solar  flare  of  28  October,  2003.  Long-lived  en¬ 
hancements  in  the  total  electron  content  (TEC)  were  observed  as  was  an  asymmetric 
increase  in  electron  temperature  due  to  increased  photoelectron  heating.  Increases  in 
electron  temperature  of  more  than  500  K  were  observed  in  the  mid-latitudes  of  the 
Northern  Hemisphere  without  a  corresponding  increase  in  the  Southern  Hemisphere. 
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Ionospheric  Response  to  Solar  Flares 
Using  an  Improved  Version  of  SAMI2 

I.  Introduction 

1 . 1  Motivation 

As  the  Air  Force  enters  the  21st  century,  its  exploitation  of  the  electromagnetic 
(E-M)  spectrum  continues  to  grow.  Long-distance  communications,  radar,  command 
and  control  of  unmanned  aerial  vehicles  (UAVs),  intelligence  gathering,  and  naviga¬ 
tion  using  the  Global  Positioning  System  (GPS)  must  all  pass  through  or  reflect  off 
the  region  of  the  atmosphere  called  the  ionosphere.  Changes  to  the  structure  of  the 
ionosphere  due  to  solar  flares  can  alter  or  prevent  the  propagation  of  E-M  waves  caus¬ 
ing  errors  in  single-frequency  GPS,  over-the-horizon  radar,  and  High  Frequency  (HF) 
radio  communications  from  D  and  E  region  absorption  of  radio  waves.  In  order  to 
anticipate  and  mitigate  the  effects  of  solar  flares  on  Air  Force  systems,  the  response 
of  the  ionosphere  must  first  be  modeled  and  characterized. 

During  a  solar  flare,  the  solar  output  shortward  of  200  nm  can  vary  by  several 
orders  of  magnitude,  on  time  scales  ranging  from  seconds  to  tens  of  minutes.  This 
energy  is  preferentially  absorbed  at  low  altitude  in  the  atmosphere,  between  60  and 
300  km,  increasing  electron  (and  ion)  density  and  electron  temperature,  which  are 
further  modified  by  the  neutral  winds,  electric  fields,  and  the  Earth’s  magnetic  field 
[Schunk  and  Nagy,  2000]. 
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1.2  Problem  Statement 


Most  ionospheric  models  are  not  adequate  for  the  study  of  solar  flares,  which 
occur  on  time-scales  of  minutes  to  hours.  The  main  limitation  is  insufficient  temporal 
resolution  of  the  solar  irradiance.  Many  current  models  use  the  Solar  Extreme  Ultravi¬ 
olet  Flux  Model  for  Aeronomic  Calculations  (EUVAC)  [ Richards  et  al,  1994]  to  specify 
the  solar  irradiance.  This  model  uses  the  10.7  cm  solar  flux  (10~22W  m-2  Hz^1),  a 
daily  measurement,  as  a  proxy  for  the  solar  irradiance.  This  yields  a  solar  spectrum 
which  is  constant  over  the  course  of  a  day,  only  varying  as  a  function  of  solar  zenith 
angle  (y)  and  optical  depth  ( tau ).  A  solar  irradiance  spectrum  which  varies  on  the 
order  of  60  seconds,  or  less,  is  necessary  to  capture  the  detail  of  a  solar  flare,  which 
begins,  peaks,  and  decays  over  a  period  of  a  few  minutes  to  an  hour.  Over  the  past 
decade,  new  satellite  measurements  of  solar  irradiance  have  allowed  the  development 
of  a  more  accurate  model  of  solar  flare  irradiance  at  a  cadence  sufficient  for  modeling 
the  ionospheric  response  to  a  flare.  Additionally,  the  processing  power  of  personal 
computers  has  advanced  to  the  point  where  running  a  research-quality  ionospheric 
model  no  longer  requires  a  large  investment  in  hardware. 

This  research  takes  an  ionospheric  model  from  the  Naval  Research  Laboratory 
(NRL),  Sami2  is  Another  Model  of  the  Ionosphere  (SAMI2)  [Huba  et  al,  2000],  and 
makes  several  modifications  to  it.  First,  the  model  code  is  rewritten  to  use  the  22-bin 
solar  irradiance  scheme  from  Solomon  and  Qian  [2005].  This  allows  the  inclusion, 
by  the  use  of  multiplicative  factors,  of  additional  ionization  clue  to  dissociative  and 
photoclectron  ionization.  Next,  the  photoelectron  heating  parameterization  is  re- 
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placed  with  a  new  one,  which  responds  to  the  increased  short-wavelength  irradiance 
during  a  flare.  Finally,  further  modification  enables  the  model  to  move  from  using 
EUVAC  irradiance  to  the  high-cadence  irradiance  of  the  Flare  Irradiance  Spectral 
Model  (FISM)  [ Chamberlin  et  al,  2007].  Once  these  changes  are  in  place  the  model 
is  used  to  investigate  the  ionospheric  response  to  a  flare. 

1.3  Previous  Work 

Studies  of  the  ionospheric  response  to  solar  flares  have  been  conducted  by  several 
researchers,  however,  most  focus  on  the  largest  (and  rarest)  events.  This  bias  is 
perfectly  reasonable  since  these  large  events  provide  the  least  ambiguous  effects,  and 
this  approach  will  be  taken  here.  The  ionosphere  occupies  a  region  of  the  atmosphere 
that  is  difficult  to  measure,  either  in-situ  or  remotely.  It  is  located  too  high  for 
balloon-borne  instruments  and  the  atmospheric  density  is  too  high  for  spacecraft  to 
survive  for  long,  especially  below  200  km.  Although  dual-frequency  GPS  can  be 
used  to  measure  electron  density,  this  is  an  integrated  measurement  which  yields  no 
information  about  variations  with  height.  On  the  other  hand,  both  incoherent  scatter 
radar  and  the  ionospheric  sounder  (ionosonde)  can  probe  the  ionosphere  from  the 
ground,  yielding  information  about  electron  density  and  temperature;  however,  the 
former  depends  upon  concurrence  of  data  campaigns  and  flares  ,  and  the  latter  is  often 
hampered  by  D  and  E  region  absorption  during  a  strong  flare.  The  consequence  is 
that  relatively  few  examples  exist  of  the  modeled  response  to  a  flare  being  compared 
to  real  data. 
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Two  researchers  have  compared  the  ionospheric  response  to  flares  using  ionosonde 
data.  Smithtro  et  al.  [2005]  and  Parsons  [2006]  both  investigated  the  ionospheric  re¬ 
sponse  to  flares  with  the  peak  flux  from  10-5  — 10~4  W/m2  (GOES  M— XI  class) 
using  the  1-D  Time-Dependent  Ionospheric  Model  (TDIM),  developed  at  Utah  State 
University.  In  order  to  specify  the  flare  irradiance  Parsons  [2006]  used  FISM  [ Cham¬ 
berlin  et  al,  2007]  in  place  of  the  EUVAC  [Richards  et  al. ,  1994],  The  Global  Airglow 
(GLOW)  model  [ Parsons ,  2006]  was  used  to  calculate  photoelectron  ionization  and 
heating  rates.  Smithtro  et  al.  [2005]  used  scaling  factors  with  data  from  the  Solar 
EUV  Experiment  (SEE)  aboard  the  NASA  Thermosphere,  Ionosphere,  Mesosphere, 
Energetics,  and  Dynamics  (TIMED)  satellite  instead  of  FISM  for  the  irradiance  and 
used  a  heating  factor  based  on  the  GOES  X-ray  flux  to  drive  the  electron  volume  heat¬ 
ing  rate.  The  results  of  these  two  papers  were  contradictory,  although  not  conclusive. 
Smithtro  was  able  to  recreate,  and  pose  a  reasonable  theory  to  explain,  a  decrease 
in  the  density  at  the  F2  peak  (NmF2),  and  simultaneous  increase  in  total  electron 
content  (TEC),  just  following  the  flare  peak,  where  one  TEC  unit  (TECU  =  1016  m-2. 
Parsons  was  not  able  to  recreate  this  flare  notch,  but  showed  the  more  conventional 
increase  in  TEC  and  increase  in  NmF2. 

Several  studies  have  been  published  using  SAMI2  and  SAMI3  by  its  principal  de¬ 
signers,  Ptuba  et  al.  [2000,  2002,  2003,  2005],  including  one  on  the  effects  of  the  Bastille 
Day  solar  flare  of  July  14,  2000.  Pluba  et  al.  [2005]  used  SAMI3,  a  3-dimensional  ver¬ 
sion  of  SAMI2  which  allows  zonal  transport  of  plasma.  Instead  of  using  EUVAC,  they 
generated  the  irradiance  spectrum  via  a  time-dependent  hydrodynamic  code  [Mariska 
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et  al,  1989]  with  the  UV  and  EUV  wavelengths  limited  by  the  data  from  the  Tran¬ 
sition  Region  and  Coronal  Explorer  (TRACE)  satellite  [Huba  et  al,  2005].  Results 
were  consistent  with  previous  results  from  Tsurutani  et  al.  [2005]  and  Meier  et  al. 
[2002]  in  TEC  increases  at  and  around  the  sub-solar  point.  They  also  noted  that  at 
the  sub-solar  point,  HmF2  decreased  by  nearly  100  km  and  NmF2  increased  by  around 
10%  [ Huba  et  al,  2005].  The  authors  pointed  out  two  important  differences  in  the 
model  output  when  compared  to  the  Global  Positioning  Satellite  (GPS)  data.  First, 
the  model  predicted  a  slower  TEC  rise  time  than  was  actually  observed.  The  authors 
allow  that  there  may  be  two  causes  for  this:  SAMI3,  like  SAMI2  doesn’t  include 
ionization  by  photoelectrons,  and  the  temporal  resolution  of  the  modeled  irradiance 
(5  min)  may  be  too  slow  [ Huba  et  al,  2005].  Second,  they  found  a  persistence  in  the 
TEC  enhancement  not  reflected  in  the  GPS  data.  They  suggested  that  this  could  be 
due  to  day-to-day  spatial  variations  of  the  equatorial  anomaly  which  are  not  included 
in  the  model  [ Huba  et  al,  2005]. 

Meier  et  al.  [2002]  used  SAMI2  to  investigate  the  Bastille  Day  flare  of  July 
14,  2000,  but  replaced  EUVAC  with  a  more  suitable  irradiance  model.  A  modified 
version  of  the  Naval  Research  Laboratory  Extreme  Ultraviolet  model  (NRLELIV) 
[  Warren  et  al,  2001;  Lean  et  al,  2003]  was  used  in  place  of  EUVAC,  using  the  At¬ 
mospheric  Ultraviolet  Radiance  Integrated  Code  (AURIC)  used  to  calculate  energy 
deposition  rates,  photoelectron  fluxes,  excitation  rates,  and  column  emission  rates 
[Meier  et  al,  2002],  They  found  an  asymmetric  response  between  the  Northern  and 
Southern  Hemispheres,  with  a  40%  increase  in  NmF2  in  the  Northern  Hemisphere 
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versus  a  nearly  30%  increase  in  the  Southern  Hemisphere,  although  no  observations 
were  available  to  corroborate  the  model  results.  Their  results  were,  however,  consis¬ 
tent  with  the  work  of  Tsurutani  et  al.  [2005],  who,  for  the  same  flare,  found  that  the 
flare  subsolar  point  was  the  point  of  greatest  increase  in  TEC,  about  29  TECU  above 
the  preflare  level.  Tsuritani’s  study  focused  on  an  examination  of  the  October  28, 
2003  solar  flare,  which  was  one  a  several  strong  flares  produced  from  late  October 
through  early  November.  He  concluded  that  the  October  28  flare  produced  a  strong 
Extreme  Ultraviolet  (EUV)  enhancement  in  addition  to  the  more  classic  X-ray  peak. 
This  difference  explained  a  long-lived  (three  hours)  increase  in  the  TEC  that  was  not 
consistent  with  ionization  from  a  large  flare  clue  to  X-rays  alone  [ Tsurutani  et  al, 
2005].  The  enhancement  they  found  is  consistent  with  the  results  of  this  research 
effort,  which  is  discussed  in  greater  detail  in  Section  4.4. 
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II.  Background 


This  chapter  covers  the  principles  necessary  to  understand  the  challenges  asso¬ 
ciated  with  ionospheric  modeling  as  it  applies  to  the  the  study  of  solar  flares. 
First,  the  chapter  will  briefly  describe  the  ionosphere  and  its  structure,  followed  by 
discussions  of  ionospheric  modeling  and  the  SAMI2  model.  Next,  a  primer  on  solar 
flares  and  the  utility  of  the  F1SM  model  over  other  models  as  the  source  for  solar 
irradiance.  Finally,  an  overview  of  photoelectrons  and  electron  heating  will  cap  the 
chapter. 

2.1  The  Ionosphere 

The  ionosphere  is  the  region  of  the  Earth’s  atmosphere  consisting  of  weakly 
to  partially  ionized  plasma.  This  region  occupies  the  atmosphere  from  roughly  60  - 
1000  km;  however,  the  upper  boundary  can  range  between  600  -  1600  km  [ Schunk 
and  Nagy,  2000].  A  weakly  ionized  plasma  is  one  in  which  collisions  with  neutral 
particles  dominate.  Ion-ion  and  ion-electron  collisions  can,  therefore,  be  neglected 
in  the  momentum  equation.  In  a  partially  ionized  plasma,  collisions  between  neutral 
constituents  and  charged  particles  must  be  accounted  for  and  in  a  fully  ionized  plasma, 
such  as  that  of  the  protonosphere,  the  neutral  density  is  low  enough  that  collisions 
between  neutral  atoms  and  charged  particles  can  be  ignored  [ Schunk  and  Nagy,  2000]. 
The  terrestrial  ionosphere  is  divided  into  four  regions;  in  order  of  ascending  altitude 
they  are  the  D,  E,  Fi,  and  F2  regions.  The  different  regions  are  identified  by  relative 
maxima  and  are  primarily  distinguishable  by  electron  density,  production  and  loss 
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mechanisms,  and  persistence  in  response  to  diurnal  changes  in  solar  ionizing  radiation. 
The  D  region,  existing  from  60  -  100  km,  is  composed  largely  of  neutral  gas  with 
a  small  ionized  component  of  molecular  ions.  Therefore,  it  is  weakly  ionized  and 
collisions  with  ions  can  be  ignored.  Due  to  its  high  neutral  density,  chemical  processes 
dominate  the  evolution  of  plasma  in  this  region.  When  photochemistry  dominates, 
the  transport  term  in  the  ion  continuity  equation  can  be  dropped  yielding: 


dns 

dt 


Ps  ~~  Ls 


(2.1) 


where  Ps  and  Ls  are  production  and  loss  mechanisms  for  each  ion  species  and  ns  is  the 
density  in  cm-3.  Production  terms  include  photoionization,  photoelectron  ionization 
and  chemistry,  while  the  main  loss  process  is  chemistry  via  recombination.  Modeling 
the  D  region  requires  specialized  models,  owing  to  its  complicated  chemistry.  This 
region  plays  a  key  role  in  the  absorption  of  HF-radio  signals,  especially  very  low 
frequency  (VLF)  signals  used  in  military  communications;  however,  the  D  region  is 
not  treated  by  SAMI2  and  is  therefore,  not  discussed  further. 

The  E  region  is  found  between  100  -  150  km  where  the  neutral  concentration 
ranges  from  1013  cm-3  to  1011  cm-3  compared  to  ion  concentrations  on  the  order  of 
104  cm-3  to  105  cm-3.  Collisions  with  neutral  species  cannot  be  ignored;  therefore, 
like  the  D  region,  it  is  weakly  ionized.  Here,  NO+  and  Oj  are  the  dominant  ions, 
and  to  a  small  extent,  N^.  Production  and  loss  mechanisms  include  photoionization 
by  extreme  ultraviolet  (EUV)  radiation  shortward  of  102  nm,  ion-neutral  chemical 


reactions,  and  recombination.  Since  photochemistry  dominates  in  this  region  too, 
ion  transport  can  be  neglected  when  calculating  the  electron  density;  however,  it  is 
included  in  SAMI2.  Although  the  ion  concentration  in  the  E  region  generally  decreases 
at  nighttime,  it  is  often  maintained  to  some  extent  due  to  the  presence  of  long-lasting 
metallic  ions  from  meteor  ablation,  and  ionization  from  starlight.  Above  the  E  region 
it  is  no  longer  possible  to  consider  dropping  transport  from  the  momentum  equation. 

The  F]  and  F2  regions  occupy  the  atmosphere  from  150  —  250  km,  and  250  to 
~  700  km,  respectively.  The  two  regions  are  often  difficult  to  distinguish,  but  at 
night  the  Fi  region  decays  leaving  a  valley  between  the  F2  and  E  regions.  Unlike  the 
lower  ionosphere,  where  molecular  ions  are  prevalent,  in  the  F2  region  ionized  atomic 
oxygen  (0+)  is  the  most  common  ion,  followed  distantly  by  atomic  nitrogen.  The  peak 
of  the  F2  region  occurs  near  300  km;  the  region  above  300  km  is  called  the  topside 
ionosphere.  Although  1000  km  was  stated  before  as  the  upper  limit  for  the  ionosphere, 
it  is  better  defined  as  the  point  where  0+  is  overtaken  by  H+  as  the  dominant  ion. 
This  point  varies  by  latitude  and  the  level  of  solar  activity,  but  is  generally  between 
600  and  1600  km  [Schunk  and  Nagy,  2000].  The  plasma  concentration  is  roughly 
105  — 106  cm-3,  while  the  neutral  concentration  is  ~  108  cm-3.  Collisions  with  both 
neutrals  and  ions  must  be  considered,  thus  the  region  is  partially  ionized.  While  in  the 
lower  ionosphere  it  is  possible  to  neglect  transport  in  the  continuity  and  momentum 
equations,  it  is  not  possible  in  the  F2  region.  The  ion  concentration  arises  from  a 
balance  between  production,  loss,  and  transport.  The  continuity  equation,  shown 
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previously,  is  given  below  with  the  ion  transport  term  intact: 

(if)  —> 

+  V  ■  (nsus)  =  Ps-Ls  (2.2) 

where  the  production  and  loss  terms,  Ps  and  Ls ,  are  as  in  Equation  2.1,  and  V  •  ( nsus ) 
is  the  transport  term  of  the  total  derivative.  0+  is  created  primarily  through  pho¬ 
toionization  of  atomic  oxygen  by  wavelengths  shorter  than  91  nm;  however,  it  is 
mainly  lost  through  charge  exchange  reactions  with  N2  and  O2.  Unlike  the  Fi  region, 
the  F2  layer  is  maintained  at  night.  There  are  two  reasons  for  this.  First,  at  the  height 
of  the  F2  layer  the  atmospheric  density  is  too  low  for  chemistry  to  act  effectively  so 
recombination  works  slowly.  Second,  as  the  0+  density  decreases,  downward  flowing 
H+  from  the  plasmasphere  helps  maintain  the  region  through  charge  exchange  with 
neutral  O  [S chunk  and  Nagy,  2000]. 

The  plasmasphere,  or  the  protonosphere,  is  the  region  from  the  topside  iono¬ 
sphere,  out  to  several  thousand  kilometers,  that  is  dominated  by  H+  and  He+.  Like 
O,  H  is  readily  ionized  by  radiation  shorter  than  91  nm,  although  most  H+  is  actually 
generated  via  charge  exchange  with  0+.  Regardless  of  the  mechanism,  the  resulting 
ions  diffuse  upward,  constrained  to  flow  along  Earth’s  magnetic  field,  reaching  heights 
of  several  thousand  kilometers.  The  H+  ions  are  relatively  cold,  less  than  1  eV,  and 
as  a  result  the  plasmasphere  motion  is  characterized  by  co-rotation  with  the  Earth, 
rather  than  the  drift  motion  associated  with  the  higher  energy  particles  of  Earth’s 
ring  current. 
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2. 2  Ionospheric  Modeling 


Given  the  importance  of  space  weather  to  communications  and  navigation,  such 
as  D  region  absorption  of  radio  waves,  and  anomalous  propagation  due  to  variations 
in  electron  density,  it  is  not  surprising  that  ionospheric  modeling  is  an  area  of  in¬ 
tensive  investigation.  Of  particular  importance  is  the  response  of  the  ionosphere  to 
the  increased  EUV  and  X-ray  flux  from  solar  flares,  which  can  affect  the  entire  sunlit 
hemisphere.  In  order  to  realistically  model  this  response,  an  accurate  model  must  be 
mated  to  an  accurate  solar  spectrum. 

Any  model  is  an  exercise  in  balancing  computational  cost  with  accuracy.  Com¬ 
putational  cost  determines  whether  it  requires  a  120-processor  mainframe,  or  whether 
it  will  run  on  a  laptop.  Accuracy  is  not  an  easy  term  to  define,  especially  when  it 
may  take  years  of  evaluating  a  model  to  identify  its  strengths  and  weaknesses.  This  is 
true  of  meteorological  models  which  have  comparatively  easy  access  to  verifying  data, 
but  doubly  so  of  ionospheric  models  whose  domain  is  largely  inaccessible  from  direct 
measurement.  At  any  rate,  cost  and  accuracy  are  in  direct  competition.  The  most 
rigorous  models  use  the  most  comprehensive  physics  (fewest  parameterizations) ,  and 
have  a  global  (3-D)  domain.  Most  models  make  use  of  parameterizations  to  speed 
computation,  especially  in  areas  where  the  added  cost  of  a  full-physics  solution  does 
not  lead  to  significantly  better  results.  This  cost  savings  can  be  multiplied  again  if 
the  domain  is  reduced  from  three  to  two,  or  one,  dimension. 
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Accurate  solar  flare  irradiance  is  an  area  that  has  seen  recent  advancements  such 


as  FISM  [ Chamberlin ,  2005;  Chamberlin  et  al,  2007].  Irradiance  models  frequently 
use  proxies  in  order  to  specify  the  solar  output.  Rather  than  try  and  measure  the 
flare  itself,  which  is  difficult  to  measure  accurately  across  the  full  X-ray  to  EUV 
spectrum,  the  major  irradiance  models  use  some,  more  easily  measured,  index  to 
modify  a  baseline  spectrum.  These  include  the  aforementioned  FISM  and  EUVAC 
models,  as  well  as  the  Hinteregger  EUV81  (HFG)  model  [ Hinteregger  et  al.,  1981]  and 
SOLAR2000  [Tobiska,  2004]  model.  The  NRLEUV,  rather  than  being  an  empirical 
model,  is  physics-based  and  generates  the  flare  spectrum  by  modeling  a  collection  of 
coronal  flux  tubes,  each  with  its  own  temperature  and  characteristic  spectrum. 

In  order  to  properly  model  the  ionosphere,  the  following  non-inclusive  list  must 
be  considered:  direct  photoionization  and  secondary  ionization  by  photoelcctrons, 
chemistry,  plasma  dynamics,  the  neutral  atmosphere,  and  both  electric  and  magnetic 
fields.  The  continuity,  momentum,  and  temperature  equations  describe  the  evolution 
of  the  plasma,  but  can  be  simplified.  Being  a  fluid,  the  plasma  can  move,  therefore 
transport  of  must  also  be  considered  along  and  perpendicular  to  the  Earth’s  magnetic 
held. 

Models  are  frequently  grouped  into  three  types,  low-latitude,  mid-latitude,  and 
high  latitude  [ Huba  et  al.,  2000].  Since  each  regime  presents  its  own  set  of  unique 
challenges,  one  model  rarely  excels  in  all  regions.  Within  the  different  latitudinal 
regimes,  models  may  span  1,  2,  or  even  3  dimensions.  The  use  of  parameterizations  can 
ease  some  of  the  computational  cost  of  a  model  and  there  are  two  that  are  common  to 
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many  models.  The  Mass  Spectrometer  Incoherent  Scatter  (MSIS)  [ Picone  et  al. ,  2002] 
and  the  Horizontal  Wind  Model  (HWM)  [ Hedin ,  1991]  specify  the  neutral  atmospheric 
density  and  neutral  winds,  respectively.  The  use  of  these  models  greatly  shortens  the 
run  time  of  a  model  by  providing  realistic  values  without  much  overhead. 

2.3  SAMI2 

Sami2  is  Another  Model  of  the  Ionosphere  (SAMI2)  is  a  2-D  (latitude  and 
altitude)  ionospheric  model  developed  by  the  Naval  Research  Laboratory  (NRL)  [Huba 
et  al,  2000].  It  models  the  mid  and  low- latitude  ionosphere,  including  the  geomagnetic 
equator,  along  a  hemispheric  grid  from  90  —  20,  000  km.  It  was  chosen  for  this  research 
because  it  is  freely  available  for  modification  from  NRL,  and  the  2-D  nature  of  the 
model  allows  more  realistic  evolution  of  the  plasma  than  a  1-D  model. 

SAMI2  models  the  evolution  of  plasma  along  an  inter-hemispheric  magnetic  flux 
tube,  solving  the  continuity  and  momentum  equations  for  7  ion  species:  H+,  He+,  N+, 
0+,  NO+  and  N^,  and  O^.  Chemistry  is  handled  via  twenty-one  chemical  processes, 
including  recombination.  Temperatures  are  calculated  for  electrons  and  three  ions 
0+,H+,  and  He+.  The  model  also  includes  E  x  B  drift  of  the  plasma  using  either  of 
two  electric  field  models.  One  is  a  simple  sinusoidal  model,  the  other  is  from  Scherliess 
and  Fejer  [1999].  This  allows  the  plasma  to  drift  horizontally  in  the  mid-latitudes  and 
vertically  near  the  dip  equator  where  the  magnetic  held  lines  are  nearly  horizontal. 

Although  the  model  calculates  the  continuity  and  momentum  equations  directly, 
it  makes  use  of  several  important  parameterizations.  First,  EUVAC  is  used  to  specify 
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the  daytime  solar  photon  flux.  EUVAC  uses  37  wavelength  bins  based  on  a  reference 
flux  scaled  by  the  daily  Fi0.r  index,  and  an  81-day  center-weighted  average  of  the 
index.  Next,  the  neutral  atmosphere  and  winds  are  calculated  using  the  MSISE-00 
and  HWM93  empirical  models.  Finally,  photoelectron  heating  is  parameterized  using 
a  method  from  Swartz  and  Nisbet  [1972], 

The  plasma  is  modeled  using  a  non-orthogonal  Eulcrian  grid,  shown  in  Fig¬ 
ure  2.1.  The  latitudinal  extent  of  the  domain  is  a  function  of  the  field  line  geometry 
and  the  height,  above  the  geomagnetic  equator  of  the  highest  field  line.  In  this  case, 
the  altitude  of  the  maximum  field  line  is  7000  km,  which  yields  a  minimum  and  max¬ 
imum  latitude  of  38°  S  and  52°  N,  respectively.  Every  fifth  field  line  is  plotted  to 
illustrate  the  geometry,  although  the  model  utilized  228  field  lines.  The  density  of  the 
field  lines,  and  gridpoints,  increases  with  decreasing  altitude  to  increase  low-altitude 
resolution.  This  can  be  adjusted  in  the  model  to  equalize  the  spacing,  but  this  was 
left  at  the  default  setting  for  this  research.  As  the  field  lines  move  farther  from  the 
equator,  they  become  more  vertical.  This  adversely  affects  the  low-altitude  spatial 
resolution  of  the  model,  and  the  only  remedy  is  to  increase  the  number  of  field  lines. 
This  increases  the  model  run  time  to  unacceptable  levels,  taking  roughly  12  hours  to 
run  a  48-hour  simulation.  The  domain,  with  contoured  data,  is  shown  in  Figure  2.2. 
Electron  density  (ne)  in  cm-3  is  plotted  for  28  October,  2003  at  1110  UT  using  the 
model  in  its  original  configuration.  The  plot  altitude  has  been  decreased  in  order  to 
show  the  low-altitude  detail.  The  model  ran  for  36  hours,  with  the  first  24  hours  used 
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Figure  2.1.  The  SAMI2  domain  for  15°  E  with  the  altitude  of  the  maximum  held 
line  (rmax)  equal  to  7000  km.  The  lowest  held  line  (rmin)  is  150  km.  Every  fifth  held 
line  is  plotted  to  clearly  show  the  geometry  of  the  grid,  which  is  comprised  of  228 
held  lines  and  301  grid  points  per  held  line. 

to  clear  any  transients  from  the  system,  then  data  is  output  every  15  minutes  for  the 
final  12  hours. 

The  plot  shows  a  density  maximum  on  each  side  of  the  geomagnetic  equator. 
This  is  consistent  with  the  equatorial  anomaly,  a  sub-tropical  electron  density  en¬ 
hancement  caused  by  E  x  B  drift  of  plasma  at  the  geomagnetic  equator.  The  plasma 
lifts  at  the  equator  and  creates  a  fountain-like  effect  as  the  as  it  spreads  out  and  de¬ 
scends  under  the  influence  of  gravity.  The  asymmetry  of  the  density  maxima,  called 
the  winter  anomaly,  is  clue  to  differences  in  the  chemistry  of  the  summer  and  winter 
hemispheres.  Both  these  effects  are  further  modified  by  the  neutral  winds,  which  can 
push  the  plasma  up  or  down  the  held  lines. 


15 


Figure  2.2.  SAMI2  electron  density  (ne)  in  cm-3  s_1  on  28  October,  2003,  1110  UT 
at  15°  E.  Solar  conditions  are  F10.7  =  270.9  and  (F10.7)  =  188.  The  domain  is  as  in 
Figure  2.1,  plotted  in  order  to  show  low-altitude  detail. 

2-4  Solar  Irradiance 


The  ionosphere  is  formed  through  absorption  of  incoming  solar  radiation  with 
wavelengths  shorter  than  200  nm.  At  these  wavelengths,  photons  possess  sufficient 
energy  to  dissociate,  excite,  or  ionize  neutral  particles.  An  example  below  shows  the 
ionization  of  diatomic  oxygen  (O2)  by  a  photon  (hv)  with  a  wavelength  shorter  than 
102.8  nm: 


O2  T  hr  — >  Oi^  T  e 


(2.3) 


where  e  is  the  free  electron,  and  hv,  is  the  energy  of  the  photon,  which  must  exceed 
the  ionization  threshold  of  the  neutral.  In  this  case  the  minimum  energy  required  to 
ionize  O2  is  12.06  eV. 
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Before  the  space  age,  this  energy,  comprising  a  small,  but  significant  portion 
of  the  solar  spectrum,  was  missed,  since  it  does  not  reach  below  the  height  of  the  D 
region.  By  the  early  1950’s,  the  extreme  ultraviolet  and  shorter  wavelength  spectra 
had  been  measured  via  rocket-borne  sensors,  and  their  importance  to  the  ionosphere 
was  firmly  established  [ Banks  and  Kockarts ,  1973].  The  arrival  of  satellite-based 
measurements  with  their  long  residence,  we  could  more  completely  measure  the  solar 
spectrum,  and  track  its  variability  over  the  solar  cycle.  Figure  2.3  illustrates  the 
variability  of  the  extreme  ultraviolet  (EUV)  down  to  10  nm,  as  well  as  the  lack 
of  variation  above  200  nm.  At  wavelengths  shorter  than  105  nm,  the  solar  flux  can 
increase  by  more  than  125%  during  active  periods,  while  the  flux  at  longer  wavelengths 
changes  by  less  than  0.1%.  While  the  vast  majority  of  the  Sun’s  output  from  visible 
light  to  radio  waves  come  from  the  5770  K  photosphere,  the  much  hotter  (105— 106  K), 
but  much  less  dense  chromosphere  and  corona,  are  responsible  for  the  EUV  and  X-ray 
portions  of  the  solar  spectrum.  During  active  periods,  the  short-wavelength  emission 
is  enhanced  by  intense  bursts  of  energy  called  solar  flares. 

2.5  Solar  Flares 

A  solar  flare  is  a  short-duration,  localized  brightening  of  the  Sun’s  chromosphere 
and/or  lower  corona.  These  explosive  bursts  of  energy  occur  on  timescales  ranging 
from  seconds  to  minutes,  making  observation  difficult  until  the  advent  of  space-borne 
instruments.  Space-based  platforms  can  persist  for  years,  even  decades,  tracking  the 
Sun  over  the  course  of  an  entire  solar  cycle,  and  maintain  nearly  constant  vigilance. 
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Solar  Irradiance 


Figure  2.3.  Variability  of  the  solar  spectrum.  The  dark  shaded  region  shows  the 
change  from  periods  of  low  solar  activity  to  active  periods. 

These  instruments  also  orbit  at  altitudes  where  absorption  is  minimal.  This  is  espe¬ 
cially  important  for  accurate  measurement  of  solar  X-rays. 

Theories  abound  as  to  the  exact  cause  of  solar  flares,  however  nearly  all  involve 
the  rapid  conversion  of  magnetic  energy  to  particle  acceleration  by  a  process  called 
magnetic  reconnection.  Aschwanden  [2004]  offers  a  comprehensive  treatment  on  the 
subject,  which  is  too  complex  to  fully  explore  here.  Figure  2.4  shows  the  enhancement 
in  ultraviolet  (UV),  and  shorter  wavelengths,  above  the  quiet  sun  level  for  both  an 
active  Sun  and  a  strong  solar  flare.  At  the  soft  X-ray  (SXR),  or  X-ray  ultraviolet 
(XUV)  wavelengths,  those  from  0.1  —  10  nm,  the  increase  can  be  up  to  several  orders 
of  magnitude,  while  in  the  UV  it  is  more  modest.  In  the  hard  X-ray  (HXR)  range, 
from  10~3  — 10-1  nm,  the  increase  can  be  five  orders  of  magnitude. 

Solar  flares  are  classified  by  their  peak  flux,  in  W/m2  from  0.1  — 0.8  nm  (SXR 
region),  as  measured  by  the  X-ray  instrument  onboard  the  GOES  class  of  weather 
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Figure  2.4.  Enhancement  of  the  solar  spectrum  during  a  large  solar  flare.  This 
demonstrates  the  order  of  magnitude  changes  in  the  X-ray  range  while  the  visible  and 
longer  wavelengths  remain  unchanged.  Adapted  from  [Foukal,  2004], 


satellites.  Table  2.1  shows  the  current  system  used  to  classify  solar  flares,  based  on 
the  peak  flux.  Photons  in  this  energy  range  are  caused  by  electrons  accelerated  to 
relativistic  speeds  (0.20.5  c)  in  the  solar  corona  [. Aschwanden ,  2004],  These  relativistic 
particles  then  loss  their  energy  via  bremsstrahlung ,  or  braking,  collisions  with  atomic 
nuclei  in  the  transition  region  and  chromosphere.  The  plot  in  Figure  2.5  shows  the 
temporal  evolution  of  a  flare  with  idealized  precursor,  impulsive,  and  gradual  phases. 
The  precursor  stage  occurs  just  a  few  minutes  before  the  impulsive  phase,  from  0  - 
4  min,  and  is  mainly  in  the  SXR  and  XUV  ranges.  The  impulsive  phase  begins  and 
ends  within  just  a  minute  or  two  and  is  measurable  in  all  but  the  Ha  and  SXR  where 
it  is  masked  by  the  gradual  phase.  This  phase  is  also  seen  in  the  microwave  region 
and  in  the  100  -  300  MHz  as  a  Type  III  radio  burst.  The  gradual  phase  shows  a 
much  longer  rise  and  fall  time  and  is  not  easily  seen  in  the  HXR  band,  or  in  the  radio 
frequencies.  The  bottom  trace  shows  the  T  +  25  min  arrival,  at  1  AU,  of  energetic 
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protons  and  electrons.  Our  understanding  of  solar  flare  evolution  has  been  pieced 


together  over  the  past  two  decades  and  has  benefited  greatly  from  the  additional  data 
provided  by  satellite  observations.  A  number  of  satellites  are  capable  of  imaging  solar 

Tabic  2.1.  Solar  flare  classification  based  on  measured  flux  in  W  m~2  between  0.1  and 
0.8  nrn  from  the  GOES  satellite.  An  XI. 7  flare  =  1.7  x  10-4  W  nr-2. 


Classification 

Peak  Flux  0.1  -  0.8  nm  (W  m  2) 

X 

1  x  10"4 

M 

1  x  10"5 

C 

1  x  10~6 

B 

i  x  nr7 

A 

1  x  10~8 

flares  at  frequencies  ranging  in  X-ray,  and  longer,  wavelengths.  The  (Rueven)  Ra- 
maty  High  Energy  Solar  Spectroscopic  Imager  (RHESSI),  and  the  Yohkoh  Soft  X-ray 
Telescope  (SXT)  are  used  to  measure  X-rays,  while  the  Solar  and  Heliospheric  Ob¬ 
servatory  (SOHO),  and  Transition  Region  And  Coronal  Explorer  (TRACE)  satellites 
are  both  optimized  for  the  EUV  range.  The  availability  of  new  measurements  of  solar 
irradiance,  and  the  observation  at  many  wavelengths  of  solar  flares,  has  enabled  the 
creation  of  F1SM  model. 
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Figure  2.5.  Temporal  evolution  of  a  solar  flare  showing  the  precursor,  impulsive 
and  gradual  phases  in  various  emission  ranges.  Adapted  from  Foukal  [2004], 
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2.6  FISM 


The  Flare  Irradiance  Spectral  Model  (FISM),  developed  by  Chamberlin  [2005], 
is  an  empirical  model  of  solar  irradiance  that  produces  a  spectrum  from  X-ray  to 
vacuum  ultraviolet  (VUV)  wavelengths  (0.5  — 192.5  nm)  at  1  nm  resolution,  with 
a  temporal  resolution  of  60  seconds  [ Chamberlin  et  al,  2007].  Data  from  several 
platforms  were  used  to  construct  a  reference  spectrum  that  is  then  increased  based 
on  daily,  and  higher  cadence  proxies,  to  model  both  flare,  and  daily,  changes  in  solar 
irradiance.  The  spectrum  is  created  using  up  to  six  different  daily  proxies,  and  one 
solar  flare  proxy.  Some  of  the  other  solar  flux  models  currently  favored  by  researchers 
are  EUVAC,  the  Hinteregger  EUV81  (HFG)  model  [ Hinteregger  et  al.,  1981],  the 
SOLAR2000  model  [ Tobiska ,  2004],  and  the  NRLEUV  model  [ Warren  et  al,  2001; 
Lean  et  al,  2003].  The  NRLEUV  [Warren  et  al,  2001]  model  is  physics-based,  while 
the  other  three  are  empirical  models.  Chamberlin  [2005]  demonstrated  the  inability  of 
these  models  to  accurately  reproduce  the  solar  spectrum  shortward  of  100  nm  during 
solar  minimum  and  solar  maximum  conditions,  when  compared  with  measurements 
from  the  SEE  instrument  onboard  the  TIMED  satellite  [ Chamberlin ,  2005].  There 
are  two  main  reasons  for  this:  the  proxies  used  in  the  established  models,  and  the  use 
of  older,  less  representative  data  in  their  creation.  The  other  models  rely  on  older, 
less  numerous  observations  of  the  X-ray  and  VUV  spectrum,  although  SOLAR2000  is 
now  using  the  SEE  data  and  the  NRLEUV  physically  calculates  the  spectrum.  This 
results  in  a  less  representative  reference  spectrum,  and  the  models  have  difficulty 
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capturing  irradiance  changes  due  to  solar  cycle  variations  and  flares.  FISM  is  able  to 
overcome  both  of  these  obstacles. 

In  order  to  more  accurately  characterize  the  solar  irradiance,  FISM  makes  use  of 
more  recent  data  from  a  number  of  satellites:  the  SEE  instrument  onboard  TIMED, 
the  two  instruments  on  the  Upper  Atmosphere  Research  Satellite  (UARS),  and  two 
on  the  SOlar  Radiation  and  Climate  Experiment  (SORCE)  satellite.  Combined,  these 
platforms  yielded  over  20,000  measurements  spanning  more  than  half  a  solar  cycle, 
including  portions  of  the  solar  cycle  minimum  and  maximum  [Chamberlin  et  al. ,  2007]. 
Also,  FISM  uses  six  different  proxies,  to  cover  the  range  from  0.1—121.6  nm,  and  FISM 
can  use  all,  some,  or  only  one  of  the  proxies.  The  F10.7  index  was  found  to  be  the 
least  accurate  proxy,  but  is  available  back  to  1947  and  for  the  foreseeable  future  since 
it  is  not  a  space-based  measurement. 

2.6.1  FISM  Solar  Flare  Modeling.  Although  FISM’s  default  60-second  ca¬ 
dence  make  it  ideal  for  use  in  modeling  solar  flare  irradiance,  the  accuracy  of  the 
algorithm  is  not  what  it  could  be  due  to  the  low  number,  and  variety  of,  observed 
flares  [ Chamberlin ,  2005].  Although  the  SEE  reference  data  set  contained  observations 
of  over  200  flares,  only  39  were  used  in  constructing  the  algorithm;  26  were  viewed  in 
the  impulsive  phase  and  13  in  the  gradual  phase  [ Chamberlin ,  2005].  The  3-second 
GOES  0.1  — 0.8  nm  irradiance  was  chosen  as  the  best  proxy  for  the  gradual  and  im¬ 
pulsive  phase.  Figure  2.6a  compares  the  FISM  flare  spectrum  (black)  with  both  the 
SEE  daily  average  data  (red),  and  an  average  of  16,  10-second  measurements  (blue), 
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which  were  taken  during  the  gradual  phase  of  a  GOES  X28+  flare  on  4  November, 
2003.  The  SEE  level  3  (L3)  data  is  generated  by  merging  the  daily  average  counts 
from  the  X-ray  and  EUV  channels.  This  data  has  been  calibrated  and  corrected  to 
exclude  field  of  view  errors,  data  from  the  South  Atlantic  Anomaly  and  solar  flares 
[ Chamberlin ,  2005].  The  level  3a  (L3a)  data  is  orbit  average  of  the  same  L3  data, 
but  does  include  increased  irradiance  from  flares.  FISM  and  the  level  L3a  data  are  in 
close  agreement,  and  below  20  nm,  are  almost  two  orders  of  magnitude  greater  than 
the  SEE  average  daily  flux.  This  is  more  easily  seen  in  Figure  2.6b.  Since  the  figure 
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Figure  2.6.  (a)  -  Ratio  of  FISM  and  SEE  for  the  gradual  phase  of  an  X28+  class  solar 

flare  of  11/04/03.  The  red  line  is  the  SEE  daily  average  with  the  flare  removed  from  the 
average,  the  black  is  FISM,  the  blue  line  is  an  average  spectrum  constructed  from  16,  10- 
second  observations  during  the  flare.  (b)  -  the  ratio  of  FISM  flare  spectra  to  the  SEE 
average  daily  level  3  data  in  red  and  the  SEE  level  3  flare  spectra  to  the  daily  average  in 
black.  From  Chamberlin  [2005]. 
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compares  the  ratios  of  FISM/L3  and  L3A/L3  minus  one,  a  1:1  ratio  is  plotted  as 
zero  on  the  vertical  axis.  The  L3  and  L3A  data  appear  identical  between  112  and 
130  nm,  indicating  that  even  during  a  flare  this  wavelength  range  does  not  increase, 
but  no  further  explanation  is  given  for  this  in  the  literature.  Therefore,  not  only  is 
FISM  capable  of  producing  an  irradiance  spectrum  that  varies  realistically,  in  tandem 
with  solar  cycle  variations,  but  it  also  captures  the  rapidly  evolving  spectrum  when 
enhanced  by  solar  flares.  As  the  database  of  observed  flares  is  sure  to  grow  in  the 
future,  the  solar  flare  proxy  and  the  algorithm  on  which  it  is  based  will,  no  doubt, 
become  more  realistic. 

2.1  Ionization  and  Photoelectrons 

In  the  terrestrial  ionosphere  the  three  primary  sources  of  ionization  are:  pho¬ 
toionization,  secondary  ionization  by  photoelectrons,  and  chemistry,  often  as  a  charge 
exchange  reaction.  Photoionization  results  when  an  energetic  photon  (A  <  105  nm) 
encounters  a  neutral  gas  molecule  or  atom,  and  absorption  of  the  photon  completely 
liberates  an  electron.  This  reaction  was  illustrated  in  Equation  2.3.  In  the  second,  the 
result  is  the  same,  however,  the  cause  is  an  energetic  (or  photo )  electron  with  kinetic 
energy  exceeding  the  ionization  threshold  of  the  target.  This  is  illustrated  below  for 
0+: 

0  +  ephoto  Q+  +  ephoto  +  e  (2.4) 

where  ephoto  is  the  energetic  electron  and  e  is  the  thermal  electron  liberated  from  O. 
The  photoelectron  still  retains  enough  kinetic  energy  to  ionize  another  neutral  and 
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the  thermal  electron  becomes  part  of  the  electron  gas.  The  charge  exchange  reaction 
between  N^~  and  NO  is: 

N+  +  NO  ->  N2  +  NO+  (2.5) 

A  fourth  source  is  dissociative  ionization  (DI).  DI  occurs  when  a  diatomic  neutral  is 
dissociated,  and  one  of  the  resulting  atoms  is  also  ionized.  This  is  illustrated  below 
for  N2: 

N2  +  hu  -*■  N+  +  N  +  e  (2.6) 

where  hu  is  the  ionizing  photon.  DI  increases  the  production  rate  of  atomic  ions, 
such  as  N+  and  0+,  while  reducing  the  production  rate  of  Nj  and  O^.  Richards 
et  al.  [1994]  found  that  photoelectrons  do  not  have  significant  impact  on  DI,  but  solar 
flux  does,  causing  DI  rates  to  peak  below  150  km.  The  relative  importance  of  these 
mechanisms  varies  according  to  species,  altitude,  and  input  solar  irradiance.  Solar 
flares  disproportionately  enhance  the  XUV  wavelengths,  which  generate  more  photo- 
electrons  and  cause  more  dissociative  ionization.  Of  these  processes,  photoelectron 
ionization  and  DI  are  not  treated  in  SAMI2;  however,  the  contribution  of  the  former 
is  especially  crucial  in  modeling  the  response  of  the  ionosphere  to  solar  flares. 

Since  ionization  is  electron  production,  we  can  think  in  terms  of  photoelectron 
production.  From  Schunk  and  Nagy  [2000,  Chap.  9]  we  see  that  the  altitude,  solar 
zenith  angle,  and  energy-dependent  photoelectron  production  rate  is: 


Pe  (E,  X,  z) 


(*) 

l  S 


loo  exp[-T  (A,  x,  z)]a\ (A)ps(A,  Ei)d\ 


(2.7) 
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where  the  sum  is  over  all  energy  levels  /,  and  all  individual  species  s.  The  total  ion¬ 
ization  cross  section  at  each  wavelength  is  er*(A),  and  the  branching  ratio  is  ps(X,Ei ); 
together  they  yield  a  cross  section  specific  to  each  excited  state  as  well  as  each  species. 
The  integration  is  carried  out  from  0  to  Xsl,  which  is  the  ionization  threshold  for 
species  s,  and  the  variable  Ei  in  the  branching  ratio  is  the  energy  level  of  the  final  ion 
state.  The  energy  E  of  the  photoelectron  is  the  difference  between  the  energy  of  the 
incoming  photon,  E\  and  the  energy  level,  Ei  of  the  ion  in  its  excited  state,  minus  the 
ionization  potential  (threshold)  of  the  neutral.  Therefore,  the  rate  of  photoelectron 
production  depends  on  the  energy  of  the  incoming  photon,  and  on  the  ionization  po¬ 
tential  of  the  neutral.  The  ionization  potential  is  the  minimum  energy  required  to  free 
an  electron  from  the  nucleus,  measured  in  electron  volts  (eV).  For  example,  H  and 
O  have  ionization  potentials  of  13.60  and  13.62  eV,  respectively,  and  can  be  ionized 
by  a  photon  of  A  <  91nm.  Using  H  and  O  atoms  as  our  example,  a  91  nm  photon 
has  a  little  more  than  13.62  eV  of  energy,  therefore  something  must  happen  to  the 
excess.  Ignoring  quantum  mechanical  effects,  this  is  due  to  the  much  greater  mass  of 
the  parent;  there  is  very  little  recoil  and  almost  all  the  excess  goes  into  kinetic  energy 
of  the  released  electron.  This  photoelectron  now  continues  to  deposit  energy  into  the 
ionosphere  through  two  important  mechanisms,  secondary  ionization  and  heating  of 
the  electron  gas,  although  it  also  loses  energy  through  excitation  of  the  neutral  gas. 

Secondary  ionization,  as  previously  noted,  results  from  ionization  by  an  electron 
with  kinetic  energy  that  exceeds  the  ionization  threshold  of  the  target  neutral.  The 
local  effect  is  the  same  as  if  caused  by  direct  solar  radiation,  yielding  an  electron-ion 
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pair  or  excited  state  neutral;  however  the  spatial  distribution  of  these  effects  are  not. 
The  highest  energy  photons  penetrate  deeply  into  the  atmosphere  before  being  ab¬ 
sorbed,  far  below  100  km.  As  a  result,  high  energy  photoelectrons  are  preferentially 
created  near  this  region,  and  the  high  neutral  density  traps  much  of  the  energy  here 
due  to  collisions.  Richards  and  Torr  [1988]  found  that  photoelectron  ionization  in¬ 
creased  the  production  rate  of  and  0+  by  a  factor  of  two  near  the  production 
peak.  In  examining  the  Bastille  Day  solar  flare,  Woods  et  al.  [2003]  found  that  the 
flare  caused  a  dramatic  increase  in  the  number  of  photoelectrons,  particularly  those 
generated  by  2  —  12  nm  wavelengths. 

Figure  2.7  shows  the  output  of  SAMI2,  with  and  without  photoelectron  ion¬ 
ization.  The  data  is  for  84°  W,  0°  N,  during  solar  maximum,  at  local  noon.  NO+, 
0;]“ ,  0+,  and  e  are  shown,  and  indicate  large  relative  increases  in  density  due  to  the 
addition  of  photoclectron  ionization.  increased  25  —  45%  between  100  —  150  km. 
NO+  increased  nearly  20%  at  the  density  peak,  which  also  decreased  in  altitude  by 
5  km.  0+  increased  by  roughly  40%  near  its  peak.  This  is  consistent  with  the  findings 
of  Richards  and  Torr  [1988].  The  changes  to  the  NO+  curve  are  indirectly  caused  by 
photoclectrons  since  NO+  is  generated  through  chemical  reactions  with  N+,  0+,  Oj 
and  Nj,  rather  than  through  direct  photoionization.  Therefore,  if  the  concentration 
of  these  species  increases,  so  will  the  concentration  of  NO+.  Above  190  km,  the  den¬ 
sity  decrease  slightly,  reflecting  the  dependence  of  the  NO+  chemical  production  rate 
on  the  0+  temperature,  which  itself  is  affected  by  the  increase  in  the  0+  density. 
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Figure  2.7.  Density  of  select  ions  and  electrons,  using  the  SAMI2  ionospheric  model 
with  (solid  lines),  and  without  (dashed  lines)  secondary  ionization  by  photoelectrons. 
Conditions  approximate  high  solar  activity  for  local  noon  at  0°/— 84°  Lat/Lon. 

The  second  impact  of  photoelectrons  is  in  heating  of  the  electron  gas  via  Coulomb 
collisions  [ Smithtro  et  al. ,  2005].  Also  called  superthermal  electrons,  to  differentiate 
them  from  the  cooler  (less  than  1  eV)  background  electrons,  photoelectrons  diffuse 
along  magnetic  held  lines,  cascading  downward  in  energy  at  each  encounter  with  a 
thermal  electron,  ion,  or  neutral.  In  this  way  they  are  able  to  deposit  the  original 
energy  non-locally ,  at  altitudes  high  above  their  formation  and  even  along  held  lines 
into  the  conjugate  hemisphere  [Woods  et  al. ,  2003].  This  causes  ion  temperature 
to  increase,  ultimately  increasing  the  individual  scale  heights  H,  and  temperature- 
dependent  chemical  reaction  rates.  These  in-turn,  redefine  the  ion  proportions,  and 
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subsequently,  the  total  electron  density.  The  consummate  result  of  this  additional 
energy  then,  is  an  overall  increase  in  ion  and  electron  density. 

This  chapter  introduced  the  reader  to  several  topics,  crucial  to  understanding 
ionospheric  modeling,  and  the  importance  of  solar  flares.  The  ionosphere,  a  layer  of 
plasma  that  surrounds  the  Earth,  alters  the  propagation  of  electromagnetic  energy 
that  passes  through  it,  degrading  communications  and  GPS-based  navigation.  During 
a  solar  flare,  the  electron  density  increase,  magnifying  these  effects.  More  detailed 
modeling  of  the  ionospheric  response  to  solar  flares  will  allow  the  Air  Force  to  under¬ 
stand  and  mitigate  the  impact  to  its  systems.  The  SAMI2  2-D  ionospheric  model  was 
selected  to  investigate  the  inter-hemispheric  response  to  a  flare;  however,  it  requires 
several  modifications.  The  temporal  resolution  of  the  input  solar  irradiance  from  EU- 
VAC  is  insufficient  to  capture  the  details  of  the  solar  flare.  This  will  be  replaced  by 
irradiance  from  FISM,  whose  60-second  cadence  is  sufficient  to  capture  to  fast  rise 
and  fall  time  of  the  solar  flux,  which  drives  the  model  ionization  rate.  Photoelectron, 
and  dissociative  ionization  will  also  be  added  to  SAMI2.  Photoelectron  ionization, 
driven  by  XUV  energy,  can  increase  by  orders  of  magnitude  during  a  flare.  The  next 
chapter  focuses  on  the  specific  changes  made  to  SAMI2. 
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III.  Methodology 


This  chapter  details  the  changes  to  SAMI2,  including  the  modification  to  use 
the  output  from  FISM  and  the  tests  that  were  run.  Although  SAMI2  has  the 
temporal  resolution  to  model  the  ionospheric  response  to  a  solar  flare,  it  lacks  some  key 
components  needed  to  do  it  successfully.  A  parameterization  to  account  for  secondary 
ionization  from  photoelectrons  is  added,  as  is  the  ability  to  utilize  60  second  output 
from  the  FISM  irradiance  model  in  place  of  the  once-per-day  -F10.7  index.  Finally, 
the  SAMI2  electron  gas  heating  scheme  is  replaced  with  a  new  parameterization  from 
Smithtro  and  Solomon(S&S )  [2008]. 

3.1  22-bin  Irradiance  Conversion 

The  move  from  37-bin  irradiance  to  22-bins  was  driven  by  the  desire  to  leverage 
the  parameterization  scheme  of  Solomon  and  Qian  [2005] .  This  allowed  the  inclusion 
of  dissociative  ionization  and  secondary  ionization  by  photoelectrons.  The  reduction 
in  the  number  of  bins  also  reduces  computational  time. 

Solomon  and  Qian  [2005]  developed  a  new  low-resolution  irradiance  scheme 
( S&Q )  that  provides  more  accurate  results  when  compared  to  the  37-bin  EUVAC. 
The  EUVAC  [Richards  et  al,  1994]  uses  5  nm-wide  bins,  and  single  wavelengths  cor¬ 
responding  to  solar  line  emission,  to  represent  the  solar  irradiance  from  5  —  105  nm. 
In  contrast,  Solomon  and  Qian  [2005]  use  high  resolution  where  absorption  cross  sec¬ 
tions  change  rapidly  and  coarse  resolution  where  they  change  more  slowly,  and  cover 
the  spectrum  from  0.05—105.0  nm.  The  22-bin  scheme  was  developed  by  comparing  a 
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high- resolution  solar  spectrum  to  the  wavelength-dependent  cross  sections  of  O,  O2, 
and  N2.  Through  trial  and  error,  the  number  of  bins  and  boundaries  were  adjusted 
until  they  obtained  the  best  agreement  between  the  low-resolution  and  high-resolution 
(EUVAC  or  HFG)  schemes.  Since  total  energy  deposition,  as  a  function  of  altitude, 
is  actually  the  sum  of  many  Chapman  functions  [Schunk  and  Nagy,  2000],  multiple 
overlapping  bins  are  used  where  significant  structure  occurs  in  the  solar  spectrum. 
This  allows  a  more  accurate  representation  of  the  distribution  of  energy  with  altitude 
Solomon  and  Qian  [2005]. 

The  S&Q  method  includes  an  energy  dependent  ratio  of  secondary-to-direct 
ionization  (pe/pi)  to  parameterize  ionization  by  photoelectrons  for  each  species.  The 
ratios  were  calculated  for  each  bin,  and  species,  as  function  of  optical  depth  (r).  They 
ran  the  GLOW  model,  for  a  wide  variety  of  solar  conditions,  independently,  for  each 
flux  bin.  They  concluded  that  pe/pi  varied  little  with  solar  activity,  and  only  at  ex¬ 
tremely  high  solar  zenith  angles.  They  also  found  that  pe/Pi  was  nearly  constant  with 
optical  depth,  except  for  r~0  [Solomon  and  Qian,  2005].  Therefore,  only  one  value, 
taken  at  r  =  1  sufficed  for  each  bin  and  all  altitudes.  It  should  be  emphasized  that  this 
scheme  only  yields  the  fraction  of  additional  photoelectron  ionization  locally.  Since 
it  does  not  explicitly  calculate  fluxes  of  photoelectrons,  non-local  photoelectron  ef¬ 
fects  are  absent.  Photoelectrons  are  constrained  to  travel  along  the  Earth’s  magnetic 
field  lines,  and,  in  the  absence  of  collisions,  travel  from  one  hemisphere  to  the  other 
[Woods  et  al,  2003].  Richards  and  Torr  [1988]  studied  secondary  ionization  to  EUV 
photoionization  ratios  and  demonstrated  that  such  transport  effects  were  not  signih- 
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cant  at  altitudes  near  the  ionization  peak.  They  did  find,  however,  that  photoelectron 
flux  could  vary  by  up  to  a  factor  of  two  at  high  altitudes,  depending  on  the  season, 
due  to  the  presence  or  absence  of  electron  flux  from  the  conjugate  hemisphere. 

The  original  EUVAC  uses  a  reference  solar  flux  which  is  based  on  the  F74113 
rocket  flight  (April  23,  1974)  [ Richards  et  al,  1994],  The  reference  flux  is  scaled  for 
solar  activity  in  each  wavelength  bin  by: 

F*  =  F74113*  [1  +  Ai(P  —  80)]  (3.1) 

where  F74113*  is  the  reference  flux,  A*  is  a  scaling  factor  and  P  is  the  proxy,  given 
by: 

P  —  ((■  Fio.r)  +  Fio.t)/2  (3.2) 

Values  appropriate  to  the  22-bin  scheme  for  A,  and  the  reference  flux  are  given  by 
Solomon  and  Qian  [2005],  as  are  cross  sections  for  O,  N2,  and  02.  Cross  sections  for 
He  and  N  were  calculated  by  Smithtro  [private  communication]. 

Figure  3.1  compares  SAMI2,  using  the  old  37-bin  EUVAC  and  the  new  S&Q 
22-bin  EUVAC,  without  photoelectrons.  The  figure  shows  the  photoion  production 
rate  (cm-3  s^1)  near  local  noon  at  0.0°  N,  84.0°  W,  which  has  a  magnetic  latitude 
of  11°  N.  Solid  lines  represent  the  37-bin  model  and  dashed  lines  represent  the  new 
low-resolution  bin  structure.  The  major  differences  are  significant  increases  in  the 
production  rates  of  0+  and  Oj  below  180  km  and  125  km,  respectively.  Near  120  km, 
0+  production  roughly  doubles,  and  near  100  km  the  increase  is  by  more  than  300%, 
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but  it  is  only  greater  by  10—40%  in  the  vicinity  of  the  peak.  On  the  other  hand,  Oj 
production  increases  by  an  order  of  magnitude  between  100  and  125  km,  shifting  the 
altitude  of  peak  production  downward  by  about  20  km.  A  mid-latitude  comparison 
(not  shown)  at  40.0°  N,  84.0°  W  showed  similar  results.  This  shows  that  the  new 
scheme  works  equally  well  in  the  low  and  mid-latitudes. 
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Figure  3.1.  Ion  production  rate  (cm-3  s_1)  from  photoionization  using  the  37-bin 
EUVAC  (solid  lines)  and  22-bin  S&Q  (dashed  lines)  irradiance  schemes  in  SAMI2,  for 
four  major  species  (He+  not  shown).  Solar  conditions  approximate  solar  maximum 
(-^10.7  =  187),  near  local  noon  for  0.0°  N/84.00  W. 


The  increase  in  0+  production  is  primarily  due  to  the  increased  ionization  by 
wavelengths  shortward  of  5  nm,  which  are  not  included  in  EUVAC.  For  0^~,  the  main 
cause  is  increased  ionization  by  wavelengths  greater  than  91.4  nm,  where  the  O2  cross 
sections  are  larger  than  those  in  the  EUVAC  and  the  other  species  have  cross  sections 
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of  zero.  This  figure  does  not  reflect  the  importance  of  dissociative  ionization,  which 
for  A <3  nm,  would  increase  the  production  of  0+  rather  than  O^. 

It  should  be  emphasized  that  the  original  version  of  SAMI2,  using  the  37-bin 
EUVAC,  uses  ionization  cross  sections  for  O2  and  N2  which  are,  in  fact,  the  full 
absorption  cross  sections.  In  reality,  a  fraction  of  the  total  absorption  cross  sections 
result  in  dissociative  ionization,  yielding  0+  and  N+.  This  causes  the  original  version 
to  over-estimate  production  of  Oj  and  Nj,  and  underestimate  0+  and  N+.  For 
example,  the  ionization  cross  section  for  from  10  -  15  nm  is  3.81  x  10-18  cm2  in 
the  original  SAMI2;  however,  the  correct  cross  section  should  be  2.35  x  10-18  cm2, 
with  the  difference,  1.46  x  10-18  cm2,  being  the  cross  section  for  dissociative  ionization 
of  02  resulting  in  0+  [Schunk  and  Nagy,  2000].  It  was  deemed  more  prudent  to  begin 
the  22-bin  conversion,  using  the  correct  ionization  cross  sections  and  ignoring,  for  the 
present,  dissociative  ionization.  Even  using  cross-sections  that  should  over-estimate 
the  production  of  and  Nj  while  under-estimating  0+  and  N+,  the  37-bin  version 
still  under-estimates  O^"  production  below  110  km  when  compared  to  the  22-bin 
version. 

Figure  3.2  shows  the  resulting  ion  densities  (cm-3)  for  the  same  model  runs.  As 
in  Figure  3.1,  the  37-bin  and  22-bin  schemes  are  represented  by  solid  and  dashed  lines, 
respectively.  Consistent  with  the  minimal  change  observed  in  the  production  rate  of 
N+,  the  N+  density  shows  only  slight  changes,  while  the  decrease  in  production  rate 
of  N+  is  clearly  reflected  in  the  decrease  in  density.  The  electron  density  is  driven  by 
the  changes  in  0+  above  155  km  and  Oj]~  below,  the  latter,  increasing  by  a  factor  of 
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1.7  and  generating  a  clearly  defined  E  region  peak  at  107  km.  Near  and  above  the 
F2  region  peak  at  350  km  there  is  virtually  no  difference  between  the  schemes,  so  this 
region  of  the  plot  is  not  shown. 


Log10  Concentration  (cm  3) 

Figure  3.2.  Ion  density  (cm-3)  from  SAMI2  using  the  37-bin  EUVAC  (solid  lines) 
and  22-bin  S&Q  (dashed  lines)  irradiance  schemes.  Five  major  ion  species  are 
shown.  Solar  conditions  approximate  solar  maximum  (Fj0.7=  187),  near  local  noon 
for  0.0°  N/84.0°  W. 
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3. 2  Addition  of  Dissociative  and  Photoelectron  Ionization 

Two  vital  components  of  the  S&Q  scheme  are  the  parameterizations  that  ac¬ 
count  for  secondary  ionization  by  photoelectrons  and  dissociative  ionization  (DI). 
Enabling  SAMI2  to  take  advantage  of  these  required  reworking  the  portion  of  the 
model  code  that  calculates  the  photoionization  rate  for  each  species.  In  the  original 
version,  the  photoionization  rate  q,  of  species  /,  at  altitude  z,  is  calculated  by: 


qfz)  =ni(z)^2<jl(  A)  Fj(A)  exp 

A 


a 


oam{ A)  /  nm(s)ds 


b 


(3.3) 


where,  nfz)  is  the  density  of  the  neutral  species  l  at  z,  and  a\  is  the  ionization  cross 
section  of  the  Ith  neutral  [ Huba  et  al. ,  2000].  Part  b  yields  the  ionizing  flux  at  height 
z,  based  on  the  incident  solar  flux,  F*(A),  and  the  absorption  of  the  flux  based  on 
the  absorption  cross  section  cr^  of  species  m.  Only  O,  02,  and  N2  are  treated  as 
significant  absorbers  of  the  incoming  flux. 

In  the  new  calculation,  part  b  is  essentially  unchanged,  however,  the  calculation 
is  now  held  in  an  array  of  the  flux  for  each  wavelength  bin  at  height  z  shown  below: 


Qi(z)  =  ^2ni(z)<r}{X)(t>g{x)  (3-4) 

A 

where  4>Z(X)  irradiance  at  height  z.  Moving  nfz)  inside  the  summation  allows  straight¬ 
forward  calculation  of  photoelectron  ionization  and  DI  by  applying  of  the  appropriate 
ratios  from  Solomon  and  Qian  [2005]  to  the  calculated  photoionization  rate.  This 
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yields: 

q?ta\z)  =  tf(z)  +  q*(z)  +  (*)  (3.5) 

where  q f*  is  direct  photoionization,  is  ionization  rate  due  to  photoelcctrons,  and 
qfDI(z)  is  the  rate  due  to  DI.  Multiplying  Equation  3.4  by  the  secondary  ionization 
ratio  Pe/'Pii  yields  the  ionization  rate  due  to  secondary  ionization,  while  DI  is  calcu¬ 
lated  by  using  the  appropriate  cross  section  for  a\.  Below  is  an  example  of  calculating 
the  secondary  ionization  rate: 

9?‘(*)  =  ^2Pe/PiWni(z)al(X)(/)z(X).  (3.6) 

A 

Before  adding  photoclectron  ionization,  it  is  appropriate  to  first  discuss  the 
impact  of  properly  treating  DI.  Figure  3.3  shows  the  ion  production  rate  (cm”3  s”1) 
for  the  S&Q  scheme.  The  solid  lines  indicate  the  production  rate  using  the  absorption 
cross  sections  (as  in  the  original  37-bin  SAMI2)  for  O2  and  N2,  instead  of  the  smaller 
O2  and  N^"  ionization  cross  sections  given  by  Solomon  and  Qian  [2005].  0+  and 
N+  are  calculated  using  the  same  ionization  cross  sections  as  in  Section  3.1,  however 
there  are  no  additional  ions  produced  from  the  dissociative  ionization  of  O2  and  N2, 
respectively.  The  dashed  lines  indicate  the  model  run  with  the  DI  fractions  from  O2 
and  N2  added  to  the  production  of  0+  and  N+.  At  this  point  the  corresponding 
fractions  are  now  absent  from  O^"  and  N^. 

As  expected,  properly  accounting  for  dissociative  ionization  reduces  the  pro¬ 
duction  rates  of  Oj  and  and  increases  the  rates  of  0+  and  N+.  Neglecting  DI 
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overestimates  the  production  of  Oj  and  by  factors  of  7  and  2,  respectively.  Both 
0+  and  N+  show  increases  at  the  same  altitude  as  the  greatest  decreases  for  Oj  and 
Nj.  In  the  case  of  0+,  this  results  in  a  relatively  large  increase  far  below  the  produc¬ 
tion  peak  and  only  a  slight  change  at  the  peak.  On  the  other  hand,  the  increase  for  N+ 
is  by  more  than  three  orders  of  magnitude,  lowering  the  altitude  of  peak  production 
to  139  km.  Near  the  altitude  of  the  N+  density  peak  (See  Figure  3.4)  at  300  km,  the 
increase  is  a  more  modest  factor  of  2.9. 


Figure  3.3.  Ion  production  rate  (cm-3  s_1)  for  SAMI2  using  the  22-bin  S&Q  ir- 
radiance  scheme  with  (dashed  lines)  and  without  (solid  lines)  dissociative  ionization 
for  four  major  species.  Photoelectron  ionization  is  not  included.  Solar  conditions 
approximate  solar  maximum  (^10.7  =  187),  near  local  noon  for  0.0°  N/84.00  W. 


Figure  3.4  shows  the  resulting  ion  and  electron  densities  (cm  3)  for  the  S&Q 
scheme  with  (dashed),  and  without  (solid)  the  correction  for  dissociative  ionization. 
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The  effect  on  the  electron  density  is  minimal  above  200  km.  The  decrease  in  total 
electron  density  from  ~125  to  200  km,  despite  an  increase  in  the  density  of  Oj,  is  a 
result  of  the  off-setting  decrease  in  the  density  of  NO+  (not  shown),  which  is  generated 
solely  through  chemical  reactions.  In  this  range,  the  NO+  concentration  is  an  order 
of  magnitude  greater  than  that  of  Oj,  the  decrease  of  the  former  overwhelming  the 
increase  from  the  latter.  By  185  km,  the  0+  density  is  an  order  of  magnitude  greater 
than  NO+  and  completely  dominates  the  plasma  by  200  km. 


Figure  3.4.  Ion  density  (cm-3)  for  22-bin  S&Q  without  (solid  lines)  and  with 
(dashed  lines)  proper  accounting  for  dissociative  ionization.  Solar  conditions  approx¬ 
imate  solar  maximum  (F10.7  =  187),  near  local  noon  for  0.0°  N/84.00  W. 

Figure  3.5  shows  the  result  of  including  the  parameterization  for  photoelectron 
impact  ionization  on  the  ion  production  rate  (cm-3  s-1).  The  most  dramatic  effects 
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are  seen  between  100  and  150  km,  especially  in  the  greater  production  of  Nj,  which 
increases  by  an  order  of  magnitude  from  100  — 120  km  and  by  ~  60  — 100%  in  the 
vicinity  of  the  production  peak  at  140  km.  Similarly,  0+  increases  by  an  order  of 
magnitude  at  100  km,  decreasing  to  a  100%  increase  at  100  km  and  to  a  40%  increase 
at  240  km,  which  remains  roughly  constant  up  to  400  km.  The  change  in  is  less 
pronounced,  reaching  a  peak  increase  of  40%  at  110  km.  N+  experiences  an  order 
of  magnitude  increase  from  100  to  148  km,  resulting  in  a  double-peak  profile,  with  a 
secondary  maximum  at  110  km.  Likewise,  the  production  peak  of  N^~  drops  from  150 
to  140  km,  however  it  also  developes  a  well-defined  valley  and  secondary  peak  below 
120  km. 


Figure  3.5.  Plot  of  ion  production  rate  (cm-3  s”1)  for  22-bin  S&Q  with  (dashed 
lines)  and  without  (solid  lines)  inclusion  of  secondary  ionization  from  photoelectrons. 
Dissociative  ionization  is  also  included.  Solar  conditions  approximate  solar  maximum 
(-^10.7  =  187),  near  local  noon  for  0.0°  N/84.00  W. 
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Figure  3.6  shows  the  resulting  changes  to  the  ion  and  electron  density  (cm-3) 
for  the  same  model  runs.  The  overall  electron  density  increases  by  roughly  70%  near 
110  km.  The  difference  decreases  above  this  point,  dropping  to  a  30%  increase  by 
165  km,  before  increasing  to  40%  at  the  F2  peak.  In  order  to  gauge  the  effect  of  the 
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Figure  3.6.  Ion  density  (cm-3)  for  SAMI2  using  the  22-bin  S&Q  with  (dashed 
lines)  and  without  (solid  lines)  inclusion  of  secondary  ionization  from  photoelectrons. 
Dissociative  ionization  is  also  included.  Solar  conditions  approximate  solar  maximum 
(-^10.7  — 187)  near  local  noon  for  0.0°  N/84.00  W. 

S&Q  scheme  in  SAMI2,  the  37-bin  and  final  22-bin  models  were  both  compared  to 
a  figure  from  Solomon  and  Qian  [2005].  Although  a  quantitative  assessment  was  not 
possible,  Figure  6a  of  Solomon  and  Qian  [2005]  offers  a  qualitative  comparison  to 
another  model.  Figure  3.7a  shows  the  0<]”,  0+,  NO+,  and  electron  densities  for  the 
37-bin  (solid  lines)  and  22-bin  versions  of  SAMI2.  The  latter  includes  dissociative  and 
photoelectron  ionization.  Figure  3.7b  is  taken  from  Solomon  and  Qian  [2005]  where 
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they  used  the  NCAR  Thermosphere  Ionosphere  Electrodynamics  General  Circulation 
(TIE-GCM)  model  to  test  the  efficacy  of  the  low  resolution  scheme.  They  conducted 
three  model  runs  using:  the  standard  37-bin  EUVAC  flux,  the  22-bin  low-resolution 
scheme  for  the  EUVAC,  and  22-bins  with  measurements  taken  from  the  SEE  instru¬ 
ment  in  place  of  the  EUVAC  baseline  flux.  The  model  was  run  for  27.5°  N,  165.0°  W 
at  noon  during  solar  maximum  (C107  =  188,  and  (FW  7)  =  186).  The  most  striking 
change  between  the  37  and  22-bin  versions  of  SAMI2  is  the  development  of  an  E 
region  between  100  and  120  km.  In  this  region 


(a)  (b) 


Figure  3.7.  (a)  Ion  density  (cm-3)  for  SAMI2,  using  37-bin  EUVAC  (solid  lines),  and  22-bin 

S  &  Q  scheme  (dashed  lines),  (b)  Ion  density  resulting  from  several  irradiance  schemes  applied  to 
the  NCAR  TIE-GCM  model.  Blue:  37-bin  EUVAC;  red:  EUVAC  modified  for  22-bins;  black:  22- 
bin  scheme  using  TIMED/SEE  irradiance.  Solar  conditions  are  the  same  for  SAME  and  TIE-GSM; 
Co. 7  =  188,  and  (Co. 7)  =  186  for  local  noon  (0  UT)  at  27.5°  N/165.00  E.  Adapted  from  Solomon 
and  Qian  [2005]. 
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the  electron  density  increased  by  more  than  an  order  of  magnitude,  driven  by  a 
comparable  increase  in  the  Oj  density.  At  the  F2  peak,  the  new  scheme  yields  an 
electron  density  that  is  about  130%  larger  than  the  TIE-GCM,  and  roughly  50%  larger 
than  the  original  version  of  SAMI2.  Having  made  the  transition  to  22-bins,  adding 
DI  and  photoelectron  ionization,  the  next  step  improves  the  way  SAMI2  handles 
photoclectron  heating. 

3.3  Electron  Volume  Heating  Rate 

As  discussed  previously,  the  importance  of  photoelectrons  is  difficult  to  ignore  at 
low  altitudes  where  they  are  responsible  for  a  significant  portion  of  the  total  ionization 
rate.  A  second,  equally  important  effect  is  the  heating  of  the  thermal  (ambient) 
electron  gas,  which  in  turn  heats  the  ion  gas  through  Coulomb  collisions.  Specification 
of  the  electron  (Te)  and  ion  temperature  is  crucial  to  calculating  species  scale  heights 
(H),  and  chemical  reaction  rates.  The  electron  temperature  depends  on  the  volume 
heating  rate  (eV  cm-3  s^1),  however,  due  to  its  dependence  on  the  photoelcctron 
flux,  which  is  often  parameterized,  it  is  usually  not  calculated  explicitly  [Smithtro 
mid  Solomon,  2008].  From  Schunk  and  Nagy  [2000],  the  electron  temperature  is 
calculated  by: 


dTe 

dt 


2  1 


_ 

3  nekb .  5  ds  ds 


Qen  Qei  ~ 1“  Qphe 


(3.7) 
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where  a  is  the  heat  flux  along  the  magnetic  field  line  (s)  in  the  dipole  coordinate 
system.  The  cooling  terms  Qen  and  Qei  are  due  to  electron-neutral  and  electron-ion 
collisions,  respectively,  and  Qphe  is  the  heating  due  to  photoelectrons  [Huba  et  al, 
2000].  Following  the  work  of  Millward  et  al.  [1996],  the  ionosphere  is  divided  into  two 
regimes.  The  first  is  a  low-altitude  regime  in  which  collisions  cause  local  photoelcctron 
heating.  The  second,  a  high-altitude  regime,  is  dominated  by  transport  and  heating 
is  non-local.  In  the  local  regime,  a  parameterization  from  Swartz  and  Nisbet  [1972] 
(SN)  is  used  to  calculate  the  volume  heating  rate: 


Qphe  £ Pphion 


(3.8) 


where  the  volume  heating  rate  from  photoelectrons,  Qphe,  is  determined  by  the  pho¬ 
toionization  rate  Pphion,  and  average  heating  efficiency,  e.  The  heating  efficiency  is 
given  by: 

£  =  exp  [-f(R)]  (3.9) 

where  f(R )  is  a  fourth-order  polynomial  fit  to  a  deposition  efficiency  R  given  by: 

f(R)  =  12.75  +  6.9417?  +  1.667?2  +  0.080347?3  +  0.0019967?4  (3.10) 


R 


nQ  2  +  nNa  +  0.1no 


(3.11) 
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In  the  non-local,  high-altitude  regime,  photoclectron  production  is  nil,  therefore,  the 
volume  heating  rate  is  calculated  as  an  exponential  decreasing  function  of  the  inte¬ 
grated  electron  density  along  the  held  line: 


Qi^S)phe 


b(s 

—j — neqo  exp 
bo 


(3.12) 


where  q0  is  the  heating  rate  per  electron  (eV/s),  q0  =  Q/ne  ,  calculated  at  the 
altitude,  sy,  of  the  transition  from  local  to  non-local  heating.  The  ratio  ^  accounts 
for  expansion  of  the  flux  tube  with  %  as  the  held  strength  at  the  transition  point. 
The  integral  is  calculated  from  So  — >  s,  or  from  the  transition  point  to  the  current 
position,  s.  The  transition  points  are  determined  where  R  ~  10~3,  which  limits  the 
transition  altitude  to  approximately  to  about  325  km. 

The  SN  4th-order  fit  was  replaced  with  a  new  6th-order  fit  from  Smithtro  and 
Solomon  [2008]  ( SS ),  and  a  new  deposition  efficiency: 


R* 


In 


_ ne _ 

no  2  +  '^n2  +  nQ 


(3.13) 


where  the  total  density  of  O,  no,  is  used  instead  of  10%  of  the  total  density.  Smithtro 
mid  Solomon  [2008]  demonstrate  that  the  SN  fit  no  longer  provides  the  best  possi¬ 
ble  approximation  to  a  full-physics  treatment,  especially  in  light  of  advances  in  the 
measurement  of  the  solar  XUV  flux  and  key  cross  sections,  since  it  was  proposed  35 
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years  ago.  Also,  unlike  the  SN  model,  the  SS  scheme  is  capable  of  reproducing  the 
heating  due  to  large  solar  flares. 

The  SS  method  treats  the  local-regime,  Qphe,  using  three  broad  wavelength 
bands,  0  —  55  nm,  55  —  80  nm,  and  80  —  105  mu.  The  heating  contribution  from  the 
2nd  and  3rd  bands  are  calculated  by  equation  3.8,  while  the  0  —  55  nm  band  is  found 
by  tabulating  a  running  sum  of  the  ionization  rate  times  the  average  photon  energy 
in  each  wavelength  bin,  then  multiplying  the  sum  by  a  normalized  heating  efficiency 
per  ionization: 

i 

where  i  denotes  the  individual  wavelength  bins  [ Smithtro  and  Solomon ,  2008].  The 
normalized  heating  efficiency,  e,  is  calculated  using  equation  3.9  with  the  new  sixth- 
order  fit  and  R*  from  Smithtro  and  Solomon  [2008].  As  long  as  the  held  lines  remain 
wholly  within  one  regime,  only  possible  for  near-equatorial  lines,  the  volume  heating 
rate  is  a  smooth  curve,  invoking  only  equation  3.8.  Since  most  held  lines  do  cross  the 
transition  point,  both  equation  3.8  and  equation  3.12  must  be  reconciled  across  this 
boundary. 

The  transition  between  the  local  and  non-local  regime  presents  a  challenge  due 
to  the  formation  of  a  discontinuity  as  shown  in  Figure  3.8.  It  was  found  that  the 
discontinuity  formed  at  the  transition  point  when  using  the  conjugate  transition  point 
as  the  upper  limit  in  the  integral  from  equation  3.12.  In  the  original  SAMI2,  the 
discontinuity  is  smoothed,  however,  at  low  solar  activity  levels,  the  discontinuity  is  still 
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present.  The  alternative  is  to  extend  the  upper  limit  of  equation  3.12  to  the  footpoint 
of  the  field  line.  Although  this  incorrectly  extends  the  non-local  approximation  into 
the  collisional  regime,  the  relative  increase  above  the  local  heating  is  generally  less 
than  10%  for  high  -^10.7  values.  Figure  3.8  shows  the  electron  volume  heating  rate 
(eV / ends)  versus  altitude  along  a  field  line  for  both  low  and  high  solar  activity.  The 
southern  hemisphere  portion  of  the  field  line  is  on  the  left,  the  northern  portion  on 
the  right.  At  low  solar  activity  Figure  3.8a,  the  discontinuity  at  the  transition  point 
is  visible  using  both  the  SN  and  SS  parameterizations.  The  lower  transition  height 
using  the  SN  causes  q0  in  equation  3.12  to  start  out  much  greater  than  in  the  SS 
version,  causing  the  resulting  volume  heating  rate  to  remain  higher  at  all  altitudes. 
The  greatest  difference  between  the  SS  volume  heating  method  occurs  just  below  the 
transition  point  due  to  the  contribution  from  the  non-local  conjugate  hemisphere. 
Since  the  non-local  volume  heating  rate  goes  as  an  inverse  exponential  of  the  electron 
density,  the  contribution  rapidly  decreases  as  the  altitude  approaches  the  electron 
density  peak,  at  around  230  km.  At  low  solar  activity,  the  relative  difference  at  the 
southern  transition  height  of  311  km  is  43%,  but  rapidly  decreases  to  less  than  10%  by 
220  km,  a  pattern  that  is  mirrored  in  the  northern  hemisphere.  Above  the  transition 
points,  the  maximum  relative  difference  is  around  4%,  from  1100  —  1600  km  in  both 
hemispheres,  but  is  generally  <  2%  up  to  the  maximum  altitude  of  10,000  km.  At 
very  low  altitudes  the  new  scheme  results  in  heating  rates  that  are  tens  of  percent,  to 
an  order  of  magnitude,  greater  than  SN,  especially  below  120  km. 
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In  Figure  3.8b,  at  high  solar  activity,  a  similar  pattern  exists,  however  there  are 
some  important  differences.  The  new  algorithm  is  valid  to  a  much  greater  altitude 
than  the  SN  fit  and  the  improved  6th  order  fit  generates  higher  heating  rates  in  the 
local  region.  Although  it  is  not  discernible  in  Figure  3.8b,  the  transition  altitude  is 
279  km,  which  is  the  same  as  in  Figure  3.8a;  however,  atmospheric  expansion  under  the 
increased  solar  activity  should  raise  the  transition  height.  This  change  is  seen  using 
the  SS  method.  Using  the  new  6th  order  fit  and  R  value,  the  SS  method  supports 
a  much  higher  transition  height  of  490  and  470  km,  in  the  northern  and  southern 
hemispheres,  respectively.  As  a  result,  the  SN  volume  heating  rate  transitions  to  the 
non-local  equation  too  soon,  resulting  in  a  heating  rates  that  are  too  low  above  the 
transition  point.  The  increased  heating  rates  below  the  peak  directly  result  from  the 
incorporation  of  increased  short  wavelength  irradiance  into  the  design  of  the  new  fit. 
This  results  in  up  to  a  factor  of  four  increase  in  the  volume  heating  rate  near  100  km 
and  a  90%  increase  at  110  km.  The  two  methods  were  compared  for  a  variety  of  field 
line  locations,  times  of  year,  and  levels  of  solar  activity.  The  same  general  results  held 
for  all  conditions  so  are  not  shown.  It  was  noted,  however,  that  when  the  solar  zenith 
angle,  y,  was  greater  than  around  80°,  i.e.,  one  footpoint  was  in  shadow  while  the 
other  was  near  local  noon,  the  differcences  between  the  methods  was  greatly  reduced. 
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Heating  Rate  (eV  cm'3  s'1) 


Figure  3.8.  Electron  volume  heating  rate  along  a  field  line  using  the  parameterization  of  Swartz 
and  Nisbet  [1972]  (solid),  and  Smithtro  and  Solomon  [2008]  using  the  conjugate  transition  point 
as  the  upper  limit  for  the  non-local  regime  integral  (dotted)  and  the  conjugate  footpoint  as  the 
upper  limit  (dashed).  The  time  is  near  local  noon  on  the  vernal  equinox  to  ensure  the  same  solar 
zenith  angle  in  both  hemispheres,  (a)  Conditions  approximate  solar  minimum,  (b)  Conditions 
approximate  solar  maximum. 
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3.4  Adding  Irradiance  from  FISM 


To  complete  the  modifications  two  final  changes  were  made  to  the  code.  The 
first  involved  bracketing  the  time  of  interest  with  a  higher  data  output  frequency. 
This  allowed  the  model  to  run  as  before,  with  data  output  at  15-minute  intervals, 
except  during  the  time  of  the  flare,  when  the  output  frequency  increases  to  1-minute 
intervals. 

The  second  change  modified  SAMI2  to  read  the  60  second  irradiance  data  from 
the  FISM  model,  which  was  described  in  Section  2.6.  FISM  provides  the  irradiance 
values  in  1  nm-wide  bins  from  0.5—192.5  nm,  however,  the  Solomon  and  Qian  22-bin 
photoionization  scheme  includes  irradiance  bins  of  0.05—0.4  nm  and  0.4— 0.8  nm.  The 
FISM  files,  covering  period  from  00  —  23:59  UT,  were  combined  with  the  GOES-12 
1-minute  0.04  —  0.5  nm  and  0.1  — 0.8  nm  X-ray  flux  ( NGDC ),  using  a  data  parsing 
utility  from  Smithtro  [private  communication] . 

In  order  to  minimize  any  spikes  caused  by  mixing  EUVAC  and  FISM,  the  FISM 
irradiance  at  00  UT  was  used  in  place  of  the  EUVAC  irradiance  until  the  begin-flare 
time,  specified  in  the  SAMI2  input  file,  was  reached.  The  model  would  then  update 
the  irradiance  every  60  seconds,  using  FISM,  from  that  point  until  the  end  of  the  run. 
Some  computational  time  may  be  lost  by  the  additional  calls  to  read  the  FISM  file 
after  the  end  of  the  flare;  however,  elevated  irradiance  may  persist  for  an  hour  or  more 
beyond  the  flare  peak,  which  would  impact  the  post-flare  evolution  of  the  ionosphere. 
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This  summarizes  the  major  changes  to  SAMI2,  making  it  suitable  for  analyzing 
the  effects  of  solar  flares.  Many  smaller  changes  were  made  in  the  course  of  debug¬ 
ging  and  ensuring  the  correct  implementation  of  the  aforementioned  alterations  to  the 
model  physics,  however,  these  are  not  worth  discussing  further.  This  chapter  summa¬ 
rized  the  modifications  made  to  SAMI2  in  order  to  use  it  to  model  solar  flare.  First, 
the  37-bin  EUVAC  irradiance  model  was  modified  according  to  Solomon  and  Qian 
[2005].  This  also  incorporated  parameterizations  for  dissociative  and  photoelectron 
ionization.  Next,  the  original  photoelectron  heating  scheme  was  modified  in  favor 
of  a  new  method  from  Smithtro  and  Solomon  [2008].  Finally,  SAMI2  was  modified 
to  enable  the  use  of  high-cadence,  high- resolution  irradiance  from  FISM  [ Chamberlin 
et  al.,  2007].  SAMI2  can  now  be  used  to  model  the  ionospheric  response  to  a  flare. 
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IV.  Results 


This  chapter  describes  the  results  from  using  the  modified  version  of  SAMI2  to 
model  the  ionospheric  response  to  the  X17  flare  of  28  October,  2003.  A  brief 
discussion  of  the  flare  itself  will  be  followed  by  analyzing  the  differences  over  a  broad 
hemispheric  view,  then  narrowing  the  discussion  to  focus  on  two  specific  altitude  and 
latitude  areas.  The  results  will  focus  on  the  changes  in  ionospheric  temperatures  and 
their  relationship  to  changes  in  electron  density,  and  what  can  be  inferred  about  the 
processes  that  drive  these  changes. 

4-1  X17  Solar  Flare,  28  October,  2003 

On  28  October,  2003  an  X17.2  solar  flare  erupted  from  solar  active  region  10486, 
peaking  at  1110  UT,  as  measured  by  the  GOES  X-ray  detector.  It  was  nearly  centered 
on  the  solar  disk  and  is  currently  the  fourth  largest  flare  on  record  and  was  chosen,  first 
and  foremost,  because  it  is  a  large,  well-studied  flare.  The  extreme  increase  at  X-ray 
and  XUV,  provide  a  worst-case  scenario  to  the  model,  thus  putting  the  modifications 
through  their  paces.  Also,  having  been  the  subject  of  many  studies,  it  was  thought 
that  some  of  the  previous  work  might  provide  useful  benchmarks  for  the  results  of 
this  simulation.  Figure  4.1  shows  the  ratio  of  flare  to  pre-flare  photon  flux  for  the 
high  resolution  FISM  (solid  red),  and  the  22-bin  (solid  black).  The  ratio  is  the  flux 
at  the  flare  peak  divided  by  the  flux  at  0900  UT.  Since  FISM  is  the  primary  driver  to 
the  22-bin  flux,  they  should  be  very  similar.  This  assumption  is  correct,  except  at  the 
shortest  wavelengths  (0  —  1  nm),  where  FISM  extends  coarsely.  The  inset  shows  the 
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corresponding  GOES-12  solar  irradiance  from  the  0.1— 0.8  nm  channel  with  the  times 
chosen  for  the  ratio  depicted  by  the  vertical  dashed  lines.  0900  UT  was  chosen  as  being 
representative  of  the  pre-flare  irradiance.  The  shortest  wavelengths  increase  during 
the  flare  by  more  than  an  order  of  magnitude,  and  in  the  case  of  the  0.05  —  0.4  nm 
bin,  by  4  orders  of  magnitude,  which  is  beyond  the  upper  limit  of  the  plot. 


Figure  4.1.  Ratio  of  flare  to  pre-flare  photon  flux  (photons  cm-2s-1)  for  the  X17 
flare  of  28  October,  2003.  (Inset)  Solar  irradiance  in  (W  m-2)  from  the  GOES-12 
0.1  — 0.8  nm  X-ray  channel.  Vertical  dashed  lines  indicate  the  times  (UT)  used  in 
computing  the  ratio. 


The  version  of  SAMI2  used  in  this  simulation  included  the  changes  described 
in  previous  chapter:  secondary  ionization  by  photoelectrons,  dissociative  ionization, 
Smithtro  and  Solomon  photoelectron  heating  scheme,  and  the  Solomon  and  Qian 
( S&Q )  irradiance  scheme.  The  FISM  output,  from  [ Chamberlin  et  al,  2007],  and 
GOES-12  0.1— 0.8  nm  and  0.05  —  0.4  nm  flux  data  for  28  October,  2003  are  combined 
into  the  the  22-bin  S&Q  format  via  the  parsing  program  mentioned  in  Section  3.4.  The 
model  is  run  twice,  the  difference  between  the  runs  being  that  the  flare  run  includes 
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the  60-second  FISM  irradiance,  while  the  non-flare  run  uses  only  the  FISM  irradiance 
from  00  UT  .  Alternately,  the  model  can  be  run  using  the  EUVAC  irradiance,  however, 
this  can  cause  unnecessary  transients  in  the  output. 

The  model  grid,  being  a  hemispheric  slice,  is  defined  by  the  equatorial  (0°  mag¬ 
netic)  height  of  the  lowest  and  highest  magnetic  held  lines.  The  sub-solar  point  at 
11  UT  is  ~15.0°  E,  and  the  minimum  and  maximum  held  line  heights  over  the  mag¬ 
netic  equator  (8.0°  N  geographic)  were  90  and  7000  km,  respectively.  This  equates  to 
the  endpoints  of  the  highest  altitude  held  line  at  38.0°  S  and  52.0°  N,  although  the 
output  domain  is  somewhat  less  due  to  the  curvature  of  the  held  lines.  The  dip  angle 
of  the  held  lines  in  the  mid-latitudes  dictates  a  high  number  of  held  lines  in  order  to 
obtain  adequate  low-altitude  resolution,  however,  the  model  run  time  is  largely  dic¬ 
tated  by  the  resolution  of  the  grid.  In  order  to  keep  run  times  down,  the  number  of 
held  lines  are  kept  to  228  or  less,  and  the  number  of  grid  points  per  held  line,  to  201. 
Increasing  the  number  of  grid  points  beyond  this  made  did  not  improve  the  output, 
but  did  increase  the  run  time  by  several  hours.  The  input  parameters  are  shown  in 
Appendix  A. 

4-2  Hemispheric  Differences 

4-2.1  Height  and  Peak  Density  at  the  F2  peak,  and  TEC.  The  ionosphere 
exhibits  differences  between  the  hemispheres  on  diurnal,  seasonal,  and  solar  cycle 
scales.  Despite  being  two  months  away  from  the  winter  solstice,  the  pre-hare  ion 
density  clearly  shows  the  winter  anomaly,  characterized  by  asymmetrical  maxima  in 
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the  electron  density  (ne),  shown  in  Figure  4.2a.  This  figure  shows  the  electron  density 
in  cm-3  using  the  FISM  60-second  flare  irradiance,  but  10  minutes  before  the  flare 
peak  (1110  UT/1210  L).  The  solid  line  indicates  the  height  (km)  of  the  F2  peak 
(. HmF2 ),  calculated  with  a  FORTRAN  function  supplied  by  the  model  author  [ Huba , 
2007].  When  compared  to  the  no-flare  plot  (not  shown),  the  flare  run  shows  small 
increases  in  both  ne  and  latitudinal  extent  of  the  peaks,  but  is  otherwise  representative 
of  the  no-flare  condition  of  the  ionosphere.  As  predicted  by  the  seasonal  anomaly,  in 
the  winter  (northern)  hemisphere,  the  HmF2  is  lower  and  the  density  at  the  F2  peak  is 
higher.  Also,  while  the  S.  Hemisphere  maximum  covers  a  broad  extent,  both  in  terms 
of  height  and  latitude,  the  N.  Hemisphere  maximum  is  much  more  limited.  Worth 
noting  is  the  deep  trough  in  ne  at  37°  N;  the  Log10  =  6.0  contour  more  than  300  km 
lower  than  in  the  S.  Hemisphere.  Figure  4.2  b  shows  the  electron  density  29  minutes 
after  the  flare  peak;  the  two  maxima  have  increased  substantially  in  magnitude  and, 
south  of  20°  N,  HmF2  has  lowered  at  nearly  all  latitudes.  The  odd  time  of  the  output 
is  a  function  of  the  resumption  of  low-cadence  output  at  1139  UT.  This  lowering  of 
the  HmF2  is  particularly  noticeable  near  3.0°  S  as  the  HmF2  has  lowered  from  550  km 
to  ~350  km.  This,  however,  seems  to  be  a  function  of  the  steepness  of  the  contours, 
rather  than  a  significant  change  in  the  overall  shape. 

Figure  4.3  illustrates  the  difference  between  the  flare  and  non-flare  model  runs 
at  1139  UT  in  the  model  domain  from  36°  S  to  48°  N.  Total  electron  content  (TEC) 
(1  x  1016  m-2)  and  electron  density  (cm-3)  at  the  F2  peak  ( NmF2 )  are  shown  in 
Figure  4.3a,  while  Figure  4.3b  shows  HmF2.  Dashed  lines  indicate  the  run  without 
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Figure  4.2.  (a)  Electron  density  (Log10  cm-3)  for  28  October,  2003  at  1100  UT/1200  L  from 

modified  SAMI2  using  the  FISM  60-second  irradiance.  The  solid  line  shows  the  model  height  of  the 
F2  peak  (HmF2).  (b)  The  same  at  1139  UT,  29  minutes  after  the  flare  peak  X-ray  intensity. 


flare  irradiance,  and  the  solid  line  indicates  the  model  run  with  the  flare.  It  should  be 


noted  that  SAMI2  overestimates  TEC,  by  an  order  of  magnitude  at  some  times  during 
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the  simulation.  The  reasons  for  this  are  not  clear  at  present  and  the  model  authors  are 
investigating  Huba [private  communication].  Despite  this,  HmF2  and  NmF2  are  in  the 
range  expected,  and  the  ATEC  observed  due  the  flare  is  in  line  with  the  results  from 
others  [Tsurutani  et  al.,  2005].  In  the  top  panel,  NmF2  shows  the  greatest  change 
in  the  mid-latitudes  of  the  N.  Hemisphere,  followed  by  just  north  of  the  geographic 
Equator.  The  relative  increase  from  30°  —  45°  N  averages  ~  45%,  while  the  relative 
increase  on  the  north  side  of  the  geographic  equator,  to  about  15°,  ranges  from  20  — 
35%.  In  the  S.  Hemisphere  the  change  is  a  more  modest,  and  constant,  10  —  12%  in 
the  low  and  mid-latitudes.  While  ATEC  is  uniform  in  absolute  terms,  it  is  not  in 
relative  terms.  In  the  southern  mid-latitudes,  the  increase  of  22  TECU  equates  to 
a  12.5%  relative  increase,  while  at  the  sub-solar  point,  around  13°  S,  the  27  TECU 
amounts  to  15.3%.  In  contrast,  the  N.  Hemisphere  mid-latitudes,  with  a  baseline  of 
~49  TECU  rises  by  more  than  30%  with  the  flare. 

The  change  in  HmF2  is  least  in  the  mid-latitudes  of  both  hemispheres,  where 
it  lowers  by  10  —  15  km  in  the  S.  Hemisphere  and  by  1  —  3  km  in  the  N.  Hemisphere. 
The  height  lowers  considerably  more  near  the  equator,  however,  this  is  somewhat 
misleading.  It  is  likely  that  the  area  around  the  clip  equator  (0°  magnetic)  displays 
a  high  altitude  trough  in  electron  density  because  of  depletion  from  the  equatorial 
fountain  effect,  which  in  turn  causes  the  Appleton  anomaly.  During  the  daytime, 
an  eastward  electric  held  causes  upward,  E  x  B  drift  of  the  plasma  along  the  dip 
equator.  The  plasma  then  diffuses  along  the  magnetic  held  lines,  descending  under 
the  influence  of  gravity  to  form  the  density  maxima  on  either  side  of  the  equator. 


This  continuous  upwelling  of  plasma  depletes  the  region  around  the  dip  equator  until 
it  can  fill  again  due  to  photoionization  during  the  day.  The  strong  enhancement  of 
the  shortest  wavelengths  by  the  flare  disproportionately  increases  ionization  at  the 
lowest  altitudes  replenishing  the  depleted  region  faster  than  normal.  So  the  larger 
decrease  near  the  equator  is  certainly  enhanced,  but  is  not  solely  a  result  of  the  flare. 
This  is  borne  out  by  later  plots  (not  shown)  which  indicate  a  recovery  of  the  HmF2  to 
a  non-flare  profile  by  1300  UT. 


Figure  4.3.  Top  -  Total  electron  content  (TEC)  1  x  1016  m-2  in  black,  and  the  elec¬ 
tron  density  in  cm-3  at  the  F2  peak  ( NmF2 )  in  red,  for  28  October,  2003  at  1139  UT. 
Dashed  lines  indicate  no  flare  and  solid  lines  indicate  the  flare  run.  Bottom  -  Height 
(km)  of  the  electron  density  peak  (HmF2)  for  the  same  model  runs,  date,  and  time. 


The  evolution  of  the  flare  over  the  next  three  hours  is  shown  in  Figure  4.4a  -  c, 
where  only  TEC  and  NmF2  are  shown.  TEC  recovers  fastest  in  the  N.  Hemisphere, 
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decreasing  to  less  than  a  third  of  its  peak  increase  by  1410  UT.  In  contrast,  the  S. 
Hemisphere  decreases  by  a  third.  This  trend  continues  for  several  more  hours  and  by 
1600  UT  the  flare  TEC  in  the  S.  Hemisphere,  especially  near  the  sub-solar  point,  is 
still  15—20  TECU  above  the  non-flare  level,  while  the  N.  Hemisphere  is  only  2  TECU 
above  normal.  Not  surprising  is  the  correlation  between  the  NmF2  and  TEC,  since  the 
greatest  contribution  to  the  modeled  TEC  is  in  the  region  of  the  F2  peak.  Between 
22°  and  30°  N,  both  TEC  and  NmF2  approach  the  no-flare  values  more  rapidly  than  at 
other  latitudes.  This  is  also  due  to  the  continued  ExB  drift  and  diffusion  of  plasma 
which,  although  enhanced  by  the  flare,  is  not  the  dominant  transport  mechanism  at 
these  latitudes.  The  N.  Hemisphere  NrnF2  peak  also  exhibits  a  slow  northward  drift 
due  to  the  influence  of  the  northward,  meridional  neutral  wind.  By  1700  UT  (not 
shown)  there  is  only  a  5%  difference  between  the  flare  and  no-flare  model  runs. 
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Figure  4.4.  (a)  Total  electron  content  (TEC)  in  1  x  1016  m-2  in  black,  and  density  in  cm-3 

at  the  F2  peak  (NmF-2),  in  red  for  28  October,  2003.  Dashed  lines  indicate  no  flare,  and  solid  lines 
indicate  the  flare  run.  (b)  and  (c)  show  the  same  at  1310  and  1410  UT,  respectively. 
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4-2.2  Electron  Density  and  Temperatures.  Although  it  is  not  noticeable  at  F 
region  heights,  the  electron  density  enhancement  begins  well  before  the  flare  irradiance 
reaches  X-class  levels.  This  happens  at  lower  heights,  chiefly  below  120  km,  due  to 
the  increase  in  X-ray  and  XUV  flux  from  0.05  —  1.8  mu.  Figure  4.5  illustrates  this  by 
plotting  the  first  three  bins  of  the  S&Q  22-bin  photon  flux  (|cnr2  s^1).  The  dashed 
vertical  line  indicates  the  flare  peak  at  1110  UT.  More  than  an  hour  before  the  flare 
peak,  the  flux  in  the  0.05  —  0.4  nm  bin  flux  increases  by  an  order  of  magnitude,  to 
~5  x  105,  and  then  to  2.5  x  106  by  1030  UT.  This  increase  is  mirrored  to  a  lesser 
extent,  but  in  lock-step  with  the  0.05—0.4  nm  bin,  in  the  0.4— 0.8  nm  bin.  The  trend 
of  a  diminished  precursor  flare  continues  as  wavelengths  increase,  with  the  precursor 
becoming  less  noticeable  in  the  0.8— 1.8  nm  bin.  At  longer  wavelengths,  or  bins  4—22 
(not  shown),  only  the  main  flare  is  significant. 

The  evolution  of  the  electron  density  has  a  significant  effect  on  Te  clue  to  the 
neutral  and  ion  cooling  terms,  Qen  and  Qei,  in  equation  3.7.  Figure  4.6  shows  the 
ionosphere’s  response  to  the  flare  over  17  minutes,  from  1105  —  1122  UT.  Figure  4.6a 
shows  the  ratio  of  ne  (flare/no  flare)  and  Figure  4.6b  shows  the  difference  of  electron 
temperature  (flare  -  no  flare),  or  ATC.  At  1100  UT  (not  shown)  the  first  effects  is  an 
increase  in  the  electron  density  below  120  km.  From  100—120  km  the  enhancement  is 
around  20%,  while  below  100  km  it  varies,  from  50—80%  in  the  southern  mid-latitudes 
to  more  than  200%  elsewhere.  At  1105  UT,  shown  in  Figure  4.6,  the  changes  are  much 
more  dramatic. 
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Figure  4.5.  Photon  flux  (photons  cm-2  s^1)  for  the  three  shortest  wavelength  bins 
for  the  28  October,  2003  solar  flare.  The  FISM  model  and  GOES-12  0.05—0.4  nm  and 
0.1  — 0.8  nm  flux  are  combined  to  yield  irradiance  in  the  22-bin  Solomon  and  Qian 
[20051  format.  The  dashed  vertical  lines  indicate  the  peak  of  the  flare  as  measure  by 
GOES-12. 

The  increased  flux  prior  to  the  flare  peak  increases  the  density  at  all  altitudes 
up  to  600  km,  with  the  soft  X-rays  depositing  their  energy  at  low  altitude,  leading 
to  increases  of  200  —  800%  below  120  km.  An  +80%  maxima  begins  to  appear  from 
23°  S  to  26°  N  at  175  km,  and  an  equatorward-sloping  trough  in  the  contours  appears 
from  200  —  300  km  between  20  —  28°  N.  This  trough  coincides  with  the  position  of 
the  N.  Hemisphere  ne  maximum  and  persists  until  around  1330  UT,  at  which  point 
the  only  remaining  maxima  are  above  0°  and  40°.  The  two  ridges  at  0°  and  38°  N 
correspond  to  minima  in  the  electron  density  as  shown  in  Figure  4.2a;  the  former 
is  underneath  the  equatorial  inverted  ne  trough  while  the  latter  coincides  with  the 
N.  Hemisphere  mid-latitude  ne  trough.  In  the  middle  (1110  UT)  panel  the  southern 
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low-latitude  maximum  reaches  a  mximum  of  200%  and  has  lifted  slightly.  This  trend 
continues  in  the  bottom  panel  (1122  UT)  with  a  clearly  defined  maximum  of  80%. 
The  ridging  above  this  point  is  also  more  pronounced  with  >10%  values  now  reaching 
above  400  km.  As  this  maximum  lifts,  to  an  altitude  of  500  km  by  1410  UT,  it 
drifts  northward  and  slowly  decays,  still  being  visblc  as  a  >20%  enhancement  over 
0°  at  1510  UT.  This  maximum  appears  to  be  clue  to  XUV  production  of  0+  from  a 
combination  of  direct  photoionization  and  secondary  impact  ionization. 

The  1-D  ion  density  above  3°  S  shows  that,  between  1100  and  1110  UT,  the  0+ 
density  from  180—290  km  increases  substantially;  by  more  than  a  order  of  magnitude 
between  190  and  220  km.  The  subsequent  migration  is  caused  by  a  combination  of  a 
meridional  neutral  wind  that  acts  from  south  to  north,  and  thus  moves  the  plasma 
up  the  field  lines,  and,  as  the  field  lines  become  more  nearly  horizontal,  an  upward- 
directed  E  x  B  drift  clue  to  the  eastward  equatorial  electrojet  [S chunk  and  Nagy, 
2000].  The  plasma,  initially  formed  by  increased  production  of  0+,  is  lifted  out  of  the 
collisional  regime  to  an  altitude  where  it  can  persist,  allowing  the  parcel  to  be  easily 
identified  until  around  1700  UT. 

Figure  4.6b  shows  the  response  of  Te  to  the  flare.  As  in  Figure  4.6a,  the  re¬ 
sponse  is  seen  well  before  the  flare  peak.  Te  responds  very  quickly  to  changes  in 
photoionization;  in  some  respects  the  flare  response  is  analogous  to  that  seen  at  sun¬ 
rise,  where  the  large  increase  in  Te  is  very  rapid,  decaying  as  the  electron  density 
increases  over  1  —  3  hours.  The  increase  in  photoionization  causes  the  Qphe  term  in 
equation  3.7  to  outstrip  the  cooling  and  advection  terms.  This  effect  is  first  seen  by 
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Figure  4.6.  Ratio  of  ne  and  Te  difference  for  the  flare  and  no  flare  model  runs  for  the  28 
October,  2003  flare,  (a)  shows  the  ratio  of  electron  density,  ne,  and  (b)  shows  the  difference  in 
electron  temperature,  ATe,  between  the  flare  and  no  flare  runs. 


1102  UT  (not  shown)  as  a  broad  increase  in  Te  throughout  the  model  domain,  ele¬ 
vated  slightly  between  130  and  250  km,  as  well  as  a  maximum  centered  above  35°  N 
at  500  km.  By  1105  UT,  the  N.  Hemisphere  maximum  is  clearly  defined  as  is  a  swath 
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of  ATe>200  K  from  36°  S  to  19°  N  at  250  km.  In  the  middle  panel  (1110  UT),  the 


N.  Hemisphere  maximum  has  intensified,  while  the  previous  broad  swath  has  eroded, 
replaced  now  by  a  layer  of  ATe  >  200  K  centered  at  140  km.  What  is  happening 
makes  sense  if  the  role  of  ne  and  nn  are  considered  in  our  cooling  terms. 

The  collisional  heating  terms,  Qen  and  Qei  [Schunk  and  Nagy,  2000],  from  equa¬ 
tion  3.7  are  given  by: 


Qen  =  V  —  Um3kB  (Tn  -  Te)  (4.1) 

n 

and, 

Qei  =  E  —^3 kB  (Tj  -  Te)  (4.2) 

mi 

i 

where  the  subscript  n  and  i  indicate  all  species  of  neutrals  and  ions,  respectively.  The 
species-dependent  collision  frequencies  are: 

Uen  oc  nn  (4.3) 

and, 

Vei  OC  - E  (4.4) 

ATe  2 

where  Ai  is  the  ion  molar  mass.  Therefore,  Qen  is  proportional  to  both  ne  and  nn; 
whereas,  Qei  which  dominates  above  ~200  km,  is  proportional  to  nt  x  ne ,  or  ~  n2e 
[Schunk  and  Nagy ,  2000].  Below  200  km,  nn  is  many  orders  of  magnitude  greater 
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than  ne  and  effectively  dominates  Qen.  Also,  Tn  is  not  affected  by  an  increase  in  Te; 
therefore,  as  Te  gets  larger,  Qen  becomes  more  negative. 

The  net  result  of  this  is  that,  at  altitudes  below  200  km  the  deposition  of  SXR 
energy  results  in  a  large  increase  in  ne ,  but  only  a  modest,  ~250  K,  increase  in  Te 
due  to  the  strong  coupling  with  neutrals.  This  strong  coupling  explains  why  cooling 
is  never  observed  below  200  km  as  it  is  at  higher  altitude.  The  neutral  density 
is  unaffected  by  the  flare,  and  the  increased  SXR,  drives  both  the  higher  electron 
density  and  electron  temperature. 

Above  200  km,  Tc  responds  quickly  to  the  higher  Qphe ,  however  ne  takes  longer 
to  respond,  continuing  to  increase  past  the  flare  peak.  As  the  flare  progresses  to  the 
peak  and  beyond,  Qphe  decreases,  while  Qei  increases,  which  it  continues  to  do  as  long 
as  ne  continues  to  increase. 

At  the  flare  peak,  Qphe  is  high  across  the  domain;  however,  at  low  altitudes,  ne 
increases  more  rapidly  than  at  higher  altitudes.  Due  to  the  coupling  to  the  neutrals, 
the  initial  peak  in  Te  is  not  followed  by  a  net  cooling.  Between  200  —  300  km,  ne 
increases  rapidly  (Figure  4.6a),  and  losses  due  to  Qei  quickly  dominate  Qphe,  which  is 
still  positive  at  this  time  owing  to  the  slow  decay  of  the  flare  (still  X-class  at  12UT). 

The  maximum  centered  above  35°  N  at  550  km  follows  a  similar  cycle,  but  on  a 
longer  timescale.  Looking  back  at  Figure  4.2a,  the  maximum  ATe  occurs  in  an  area  of 
low  ne.  This  is  partially  due  to  the  new  Smithtro  and  Solomon  photoelectron  heating 
parameterization.  As  seen  in  Figure  3.8,  the  use  of  the  new  parameterization  places 
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the  transition  height  very  near  500  km  for  periods  of  high  solar  activity  (496  km  for 
the  flare).  Although  there  was  not  a  significant  difference  in  electron  volume  heating 
rates  between  similar  latitudes  above  the  F2  region,  a  strong  S.  Hemisphere  bias  was 
seen  below  200  km,  where  the  maximum  increase  was  about  200%  greater  than  that 
in  the  N.  Hemisphere.  Therefore,  the  cause  for  the  asymmetric  high  altitude  heating 
seems  not  to  be  differenct  heating  rates,  but  different  electron  densities. 

Examining  ne  and  Qvhe  along  a  field  line  show  that,  near  500  km,  the  electron 
density  in  the  S.  Hemisphere  is  more  than  an  order  of  magnitude  greater  than  in 
the  N.  Hemisphere,  but  the  volume  heating  rate  is  only  twice  as  much.  A  direct 
comparison  between  35°  S  and  35°  N  revealed  that  at  500  km,  the  heating  rate  per 
electron  (eV/s)  is  more  than  280%  greater  in  the  N.  Hemisphere.  Having  established 
the  reason  for  the  large  temperature  increase,  the  evolution  is  much  the  same  as  the 
region  over  3°  S  at  250  km;  however,  it  is  much  slower.  The  additional  photoelcctron 
heating  ( Qphe )  ceases  at  the  same  time  as  for  the  aforementioned  region,  but  ne  is 
still  very  low  and  increasing  slowly.  This  is  eveident  by  the  decrease  in  ATe.  The 
temperature  continues  to  decrease,  and  by  1154  UT  (not  shown)  a  ATC  minimum  has 
formed  at  ~450  km,  above  35°  N.  Te  continues  to  drop  and  the  area  expands  to  fill 
the  region  once  occupied  by  a  positive  ATe  by  1239  UT.  The  temperatures  do  not 
recover  full  until  after  1500  UT. 

The  lack  of  a  large  Tc  increase  below  110  km  is  due  to  the  fact  that  Te  and  Tn  are 
strongly  coupled  at  low  altitudes  and  MSIS-00,  the  neutral  atmosphere  model,  does 
not  respond  to  external  forcing  of  the  temperature.  Also,  the  high  neutral  density 


(1011  — 1012  cm  3)  insures  that  any  increase  will  be  relatively  short-lived  clue  to  the 
strong  neutral  coupling. 

4-3  Mid-Latitude  versus  Low-Latitude 

This  section  will  compare  and  contrast  the  flare  response  at  a  low-latitude, 
3°  S,  and  mid-latitude,  35°  N,  location.  These  two  location  were  chosen  because  they 
represent  latitudes  that  exhibited  the  greatest  changes  in  Te  and/or  ne  due  to  the 
flare.  Although  the  S.  Hemisphere  did  show  an  increase  in  TEC,  this  was  expected 
and  doesn’t  warrant  in-depth  analysis.  In  contrast,  at  low-latitudes,  it  is  the  southern 
hemisphere  which  exhibited  the  strongest  response  to  the  flare.  At  the  flare  peak,  the 
solar  zenith  angle  at  3°  S  was  10.8°  and  at  35°  N  it  was  ~48°.  As  in  the  previous 
section,  comparison  of  the  flare  and  no-flare  model  runs  illustrates,  not  only  how  the 
response  varies  from  low  to  mid-latitudes,  but  also  the  magnitude  of  the  changes. 

Figure  4.7  plots  the  electron  volume  heating  rate  (eV  cm-3  s_1)  at  the  flare 
peak,  for  flare  and  no-flare  runs  of  SAMI2  at  3°  S  and  35°  N.  The  solid  lines  indicate 
the  heating  rate  with  flare  irradiance  and  the  dashed  lines  indicate  the  run  without. 
The  responsiveness  of  the  Smithtro  mid  Solomon  heating  parameterization  to  the 
additional  short-wavelength  flux  is  visible  by  the  increase  of  2—3  orders  of  magnitude 
near  100  km. 

Figure  4.7a  displays  a  discontinuity  that  results  from  using  different  equations 
(3.12  and  3.14)  for  the  low  and  high  altitude  collision  regimes.  Rather  than  plotting 
the  electron  volume  heating  rate  along  a  field  line,  the  heating  rates  are  plotted  above 
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specific  points.  Therefore,  they  are  the  result  of  the  heating  rates  from  as  many  field 
lines  as  intersect  a  vertical  line  above  the  location.  Since  the  transition  point  is  at 
~  500  km,  locations  at  low  altitude  near  the  equator  sample  field  lines  which  are 
entirely  within  the  collisional  regime,  and  are  calculated  by  using  equation  3.14,  until 
the  transition  point  is  reached.  This  transition  is  apparent  by  the  sharply  decreasing 
slope  at  500  km.  Above  this  altitude  the  sampled  field  lines  reflect  a  combination  of 
local  heating  and  non-local  transport.  Using  the  original  version  of  SAMI2  revealed 
the  same  discontinuity,  although  not  at  3°  S.  The  much  lower  transition  point  (about 
280  km)  causes  the  discontinuity  to  reveal  itself  closer  to  the  geomagnetic  equator 
than  in  the  new  version. 

Figure  4.8  shows  1-D  plots  of  ion  production  rate  (cm-3  s_1)  at  3°  S  and  35°  N, 
during  the  28  October,  2003  flare.  It  bears  repeating  that  these  production  rates 
include  photoionization,  DI,  and  secondary  impact  ionization.  The  results  of  the 
simulation  using  flare  irradiance  are  shown  as  dashed  lines,  while  solid  lines  indicate 
the  model  run  without  it.  Again,  the  ionization  by  XUV  radiation  is  implied  by  the 
very  large  increase  in  ion  and  electron  production  below  200  km.  In  Figure  4.8a,  the 
increase  in  e“  production  at  100  km  is  two  orders  of  magnitude,  driven  mainly  by  the 
factor  of  39  increase  in  Nj  production.  Figure  4.8b  shows  a  similar  situation  with  Nj 
leading  production,  where  increased  production  of  0+  at  180  km  explains  the  increase 
in  electron  density  seen  in  Figure  4.6. 

Photoionization  and  photoelectron  ionization  from  the  flare  lead  to  large  differ¬ 
ences  in  electron  density  with  altitude.  Regardless  of  what  causes  the  increase,  it  is 
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Figure  4.7.  Electron  volume  heating  rate  (eV  cm-3  s_1)  at  two  different  latitudes  for  flare  and 
no-flare  runs  of  SAMI2.  The  dashed  line  is  the  no-flare  model  run  and  the  solid  line  includes  the 
flare,  (a)  shows  the  volume  heating  rate  at  3°  S,  (b)  shows  the  rate  at  35°  N  for  comparison. 


the  altitude  of  formation  that  determines  its  longevity.  Figure  4.9  shows  contour  plots 
of  the  ratio  (flare/no  flare)  of  ne  with  respect  to  altitude  and  time  at  3°  S  and  35°  N. 
The  flare  peak  is  indicated  by  arrows  above  and  below  each  plot.  At  both  locations 
the  density  begins  to  increase  at  altitudes  below  the  F2  peak  at  1000  UT.  This  is 
consistent  with  the  M-class  levels  reached  prior  to  the  flare  as  seen  in  Figures  4.1  and 
4.5.  At  the  flare  peak,  Figures  4.9a  and  b  show  the  result  of  increased  ion  production 
in  the  form  of  up  to  an  order  of  magnitude  increase  in  electron  density  near  100  km. 
The  less-enhanced  area  between  the  peaks  at  100  and  ~200  km  corresponds  to  the 
valley  between  the  E  and  F  regions,  confirmed  by  analyzing  a  1-D  plot  of  electron 
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Figure  4.8.  Ion  production  rate  in  cm-3  s_1,  for  flare  (dashed)  and  no-flare  (solid)  runs  of 
SAMI2.  The  ion  production  rate  includes  photoionization,  photoelectron  impact  ionization  and 
dissociative  ionization,  but  does  not  include  production  rate  from  chemistry.  The  rates  are  shown 
for  3°  S  in  (a),  and  35°  N  in  (b)  to  illustrate  the  differences  between  low  and  mid-latitudes. 

density  (not  shown).  As  time  moves  forward  in  Figure  4.9a,  the  200  km  density  max¬ 
imum,  caused  by  increased  0+,  diffuses  upward  as  the  neutral  wind  transports  the 
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plasma  along  the  magnetic  field.  This  results  in  a  northward  drift  until  the  plasma 
becomes  affected  by  E  x  B  drift  at  the  dip  equator.  This  was  confirmed  by  running 
the  model  without  any  neutral  winds,  and  with  winds,  but  a  reduced  electric  field. 
When  the  electric  field  was  reduced,  the  movement  of  the  low-latitude  density  max¬ 
imum  was  not  as  pronounced,  but  was  still  discernible.  Without  the  neutral  winds, 
the  winter  anomaly  fails  to  develop,  resulting  in  symmetric  density  maxima  in  both 
hemispheres.  The  density  increase  due  to  the  flare  is  almost  symmetric  with  respect 
to  the  dip  equator;  there  are  slight  differences  due  to  the  smaller  zenith  angle  in  the  S. 
Hemisphere,  and  the  differences  associated  with  the  neutral  atmosphere  model.  The 
post-flare  movement  of  the  maximum  without  wind  is  upward,  without  latitudinal 
drift.  The  effect,  of  density  increasing  in  altitude  with  time,  is  seen  as  far  northward 
as  15°  N  (geographic),  but  diminishes  rapidly  between  10°  and  15°  and  by  18°  N 
it  is  gone  entirely.  In  the  S.  Hemisphere  it  is  clearly  seen  as  far  south  as  20°  S,  or 
27°  S  magnetic.  This  result  may  also  be  caused  by  enhanced  EUV  flux  from  the 
flare,  whose  slow  decay  is  evident  in  the  increased  densities  of  20—30%  beyond  14  UT 
around  100  km. 

In  the  N.  Hemisphere  mid-latitudes,  Figure  4.9b,  the  increase  is  due  more  to 
changes  in  the  plasma  scale  height  caused  by  the  increase  in  electron  temperature.  In 
the  classical  diffusive  equilibrium  equation  [Schunk  and  Nagy ,  2000]: 


1  dni  1  1  dTp 

n,  dz  Hp  Tp  dz 


(4.5) 
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where  Hp  is  the  plasma  scale  height,  equal  to  2fcgTp  ,  and  T„  is  the  plasma  tern- 
perature.  It  can  be  shown  that  in  the  limit  of  an  isothermal  ionosphere,  the  major 
ion  and  electron  densities  decrease  exponentially  with  increasing  altitude,  at  a  rate 
determined  by  Hp.  Thus,  as  the  electron  and  ion  temperatures  increase,  so  does  the 
plasma  scale  height,  which  in  turn  increases  the  slope  of  the  density  curves  at  high  al¬ 
titude.  Finally,  looking  at  how  NmF2 ,  HmF2,  and  TEC  evolve  with  time  over  specific 
locations  may  shed  some  additional  light  on  what’s  happening. 

4-3.1  Temporal  evolution  of  NmF2,  HmF2,  and  TEC.  NmF2,  HmF2,  and 
TEC  were  examined  earlier  in  the  chapter;  however,  looking  at  how  they  change  in 
time  may  solidify  some  notions  regarding  the  flare  response  to  different  regimes,  e.g. 
collisional  vs.  transport-dominated.  Figures  4.10  and  4.11  plot  NmF2  and  HmF2 
in  panel  (a),  and  TEC  in  panel  (b)  at  3°  S  and  36°  N,  respectively.  Dashed  lines 
indicate  the  model  run  without  flare  irradiance  and  solid  lines  denote  the  run  with 
flare  irradiance.  Again  the  flare  peak  is  indicated  by  a  vertical  dashed  line.  Care 
must  be  taken  when  viewing  these  figures;  although  the  ranges  are  consistent,  the 
magnitudes  are.  This  was  done  in  order  to  fit  all  the  data  into  the  plot  area. 

The  most  obvious  difference  between  the  two  locations  is  the  decrease  in  HmF2, 
which  is  191  km  at  3°  S  versus  only  ~6  km  at  35°  N.  In  the  case  of  the  former,  the 
maximum  decrease  lags  the  flare  by  14  min,  whereas  there  is  virtually  no  lag  at  35°  N. 
This  is  a  result  of  the  very  different  pre-flare  environments.  In  the  S.  Hemisphere  the 
vertical  extent  of  density  greater  than  3.09  x  106  cm-3  extends  from  460—560  km.  In 
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Figure  4.9.  Contour  plots  of  ne  ratio  (flare/no  flare)  for  modified  SAMI2.  The  flare  peak  is 
marked  by  arrows  at  the  top  and  bottom,  (a)  shows  the  volume  heating  rate  at  3°  S,  (b)  shows  the 
rates  at  35°  N  for  comparison. 


other  words,  it  is  broad  and  flat.  Increased  ionization  of  0+  below  the  peak  causes 
the  large  drop  in  altitude  of  the  peak,  even  though  the  relative  increase  in  electron 
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density  is  only  on  the  order  of  a  few  percent.  Even  as  the  HmF2  recovers,  NmF2 
slowly  increases  as  two  things  happen:  increased  ionization  continues  due  to  the  slow 
decay  of  the  flare,  and  plasma  is  transported  along  field  lines  to  higher  altitudes.  The 
broad  vertical  extent  of  0+,  seen  in  a  1-D  plot  of  ion  density,  also  accounts  for  the 
much  greater  increase  in  TEC  seen  in  Figure  4.10b.  That  the  plasma  is  in  a  region 
where  recombination  is  slow  accounts  for  the  longevity  of  the  increase  in  both  TEC 
and  NmF2. 

In  Figure  4.11a,  NmF2  responds  almost  immediately  to  the  flare  when  compared 
to  the  slower,  but  longer-lived,  increase  seen  at  3°  S.  This  difference  is  a  result  of  the 
different  mecahnisms  involved,  and  the  altitudes  of  HmF2  at  the  two  latitudes.  At  3°  S, 
E  x  B  drift  of  the  plasma  to  heights  where  transport  dominates  and  recombination 
is  slow  accounts  for  the  longevity,  and  slow  increase,  of  NmF2.  At  35°  N,  the  increase 
is  due  to  direct  photoionization  from  the  flare,  with  transport  playing  a  much  smaller 
role.  As  such,  NmF2 ,  and  the  F2  region  as  a  whole,  in  the  N.  Hemisphere  mid-latitudes 
is  at  an  altitude  where  electron  recombination  occurs  on  a  much  faster  time  scale.  This 
is  evident  in  the  much  faster  decay  of  NmF2  at  35°  N  to  near  the  no-flare  level  when 
compared  to  3°  S.  The  dependence  of  TEC  on  NmF2 ,  shown  in  Figure  4.10b  and 
4.11b,  is  such  that  the  evolution  of  TEC  mirrors  that  of  NmF2. 
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(a) 


(b) 

Figure  4.10.  NmF2,  HmF2l  and  TEC  at  3°  S  for  modified  SAMI2  showing  evolution  with  time 
for  the  solar  flare  of  28  October,  2003.  The  flare  peak  is  marked  by  the  vertical  dashed  line  at 
1110  UT.  (a)  plots  NmF2  in  cm~3  s_1  (blue  lines)  and  HmF2  in  km  (black  lines),  (b)  plots  the 
total  electron  content  (TEC)  in  TECU  =  1016  m-2.  Solid  lines  indicate  the  flare  run,  dashed  lines 
indicate  the  no-flare  run. 
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(a) 


(b) 


Figure  4.11.  Nm F2 ,  HmF2 ,  and  TEC  at  35°  N  for  modified  SAMI2  showing  evolution  with 
time  for  the  solar  flare  of  28  October,  2003.  The  flare  peak  is  marked  by  the  vertical  dashed  line  at 
1110  UT.  (a)  plots  NmF2  in  cm~3  s_1  (blue  lines)  and  HmF2  in  km  (black  lines),  (b)  plots  the  total 
electron  content  (TEC)  in  TECU  =  1016  m2.  Solid  lines  indicate  the  flare  run,  dashed  lines  indicate 
the  no- flare  run. 
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4-4  Comparison  to  other  results 


This  section  concludes  with  a  brief  comparison  of  the  results  to  some  of  the 
previous  work  mentioned  in  Chapter  I.  In  general,  the  results  of  the  model  compared 
well  with  what  other  researchers  have  found. 

The  up  to  ~30  TECU  increase  found  in  this  study  is  comparable  to  the  maxi¬ 
mum  estimated  by  Tsurutani  et  al.  [2005].  They  measured  a  25  TECU  maximum  at 
the  subsolar  point,  based  on  TEC  from  GPS  data.  The  Libreville,  Gabon  GPS  station 
recorded  a  maximum  increase  of  ~30  TECU  on  two  of  the  six  tracked  satellites.  Lo¬ 
cated  at  0°  N/9°  E,  it  was  near  the  subsolar  point,  making  it  a  good  measurement  of 
the  actual  increase.  The  results  from  SAMI2  indicate  that,  at  least  in  the  low-latitude 
region,  that  the  long-lived  increase  in  TEC  is  due  to  the  combined  effects  of  neutral 
winds,  and  E  x  B  drift  of  plasma  from  lower  altitudes. 

This  is  validated  by  examining  the  variation  of  0+  density  with  altitude.  At 
3°  S/15°  E  the  maximum  increase  in  density  compared  to  the  no- flare  run  occurs  at 
230  km  at  the  flare  peak.  This  increase  of  nearly  100%  is  easily  followed  at  subsequent 
times.  By  1124  LIT  this  maximum  has  increased  in  altitude  and  magnitude,  exceeding 
the  magnitude  of  the  no-flare  F2  peak  at  ~600  km.  This  explains  the  cause  of  the 
abrupt  descent  of  HmF2  seen  in  Figure  4.10a, 

Figure  4.12  shows  the  evolution  of  ATe  (flare  -  no-flare)  with  time  and  altitude 
for  3°  S  and  35°  N.  Figure  4.12b  displays  the  large  ATe  due  to  the  combination  of 
photoelectron  heating  and  low  density,  which  transitions  to  a  net  cooling  as  Qphe 
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decreases  as  ne  increases.  This  is  absent  in  Figure  4.12a  due  to  the  much  higher 
electron  density,  by  an  order  of  magnitude,  at  this  altitude.  Both  latitudes  exhibit 
very  similar  low  (200—300  km)  altitude  trends,  although  the  cooling  at  3°  S  is  greater 
due  to  the  higher  density  (see  Figure  4.9)  in  a  region  where  the  cooling  rate  is  still 
proportional  to  n2.  Near  150  km,  the  important  collisions  are  with  neutrals.  Tc  is 
coupled  closely  to  the  neutral  temperature  and,  is  cooled  according  to  equation  4.1. 
This  coupling  also  dampens  the  effect  of  the  surge  of  high-energy  photons  deposited 
at  this  height.  The  weaker  coupling  at  500—600  km  at  35°  N  explains  the  much  larger 
ATe,  as  Te  is  able  to  respond  freely  to  the  influx  of  energy.  This  2-D  analysis  of  ATe 
was  also  used  by  Smithtro  and  Solomon  [2008]. 

Smithtro  and  Solomon  [2008]  also  used  the  28  October,  2003  flare  to  test  the 
new  photoelectron  heating  parameterization  using  a  1-D  model  for  40°  N.  The  re¬ 
sults  from  SAMI2  were  very  similar  with  the  following  differences  observed:  First,  the 
flare-induced  positive  ATe  at  500  km  is  around  100  K  higher;  cooling  of  the  region 
is  less  rapid,  taking  30  minutes  more  in  to  reach  ATC  =  0;  finally,  cooling  was  less 
pronounced  remaining  around  100  K  warmer  than  the  1-D  model.  The  exact  magni¬ 
tude,  and  timing,  of  the  temperature  change  is  not  as  important  as  the  fact  that  the 
change  was  observed  at  all.  As  shown  earlier,  the  electron  heating  and  cooling  rates 
at  this  altitude  vary  as  n2. 
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Figure  4.12.  Temperature  difference  (ATC  in  Kelvin)  between  flare  and  no-flare  runs 
of  the  modified  SAMI2  model.  The  scenario  is  the  X17  solar  flare  of  the  28  October, 
2003.  (a)  is  for  3°  S/15°  E  and  (b)  is  for  35°  N/15°  E,  from  10-14  UT/11-15  L. 
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V.  Conclusions  and  Suggestions  for  Future  Work 


THIS  chapter  presents  the  conclusions  drawn  from  this  research.  It  also  ex¬ 
plains  the  impact  of  the  research  to  the  United  States  Air  Force  and  makes 
suggestions  for  future  work  in  this  area. 

5.1  Summary  of  Research  and  Conclusions 

The  SAMI2  model  presents  a  unique  opportunity  to  researchers.  It  is  one  of  a 
number  of  two-dimensional  ionospheric  models,  which  are  complete  enough  to  provide 
meaningful  results,  yet  efficient  enough  to  run  on  a  laptop  computer.  It  has  also  been 
designed  to  be  easily  upgradeable  and  can  incorporate  better  models  of  variables  such 
as  the  neutral  winds  and  neutral  atmosphere  as  they  become  available. 

First  and  foremost,  the  research  was  successful.  SAMI2,  in  its  original  con¬ 
figuration,  is  not  capable  of  adequately  modeling  the  ionospheric  response  to  solar 
flares  because  of  three  key  limitations.  First,  the  input  solar  irradiance,  based  on  the 
EUVAC  model,  lacks  the  temporal  resolution  needed  to  capture  the  evolution  of  the 
flare.  Second,  secondary  ionization  by  photoelectrons,  which  are  particularly  sensi¬ 
tive  to  the  enhanced  SXR  wavelengths  from  a  flare,  is  not  included.  Finally,  SAMI2 
uses  an  older  parameterization  for  photoclectron  heating,  which  is  unresponsive  to 
the  short-wavelength  bias  of  the  flare  spectrum.  These  limitations  were  removed  with 
the  Erst  being  the  addition  of  photoelcctron  ionization. 

The  impact  of  photoelectron  ionization  in  the  vicinity  of  the  ionization  peak  for 
0+,  Oj,  and  has  been  known  for  many  years  [Richards  et  al,  1994;  Woods  et  al. , 
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2003].  At  the  production  peak,  the  effect  can  be  equally  as  important  as  direct  pho¬ 
toionization  [ Richards  et  al,  1994],  Photoelectron  ionization  was  successfully  added 
to  SAMI2,  and  demonstrated  that  the  model  could  generate  an  E  region  peak,  in  ad¬ 
dition  to  an  improved  representation  of  the  F[  peak.  This  addition  came  as  a  result 
of  incorporating  the  Solomon  and  Qian  [2005]  low-resolution  irradiance  method. 

Second,  the  original  photoelectron  heating  parameterization  was  not  capable  of 
responding  to  the  enhanced  flux  from  the  flare.  It  was  replaced  with  one  capable  of 
responding  to  the  extreme  variability  of  solar  X-ray  and  XUV  energy  associated  with 
flares.  This  new  scheme  not  only  allowed  the  model  to  generate  more  realistic  quiet 
sun  electron  heating,  but  also  to  respond  to  flare  irradiance  with  rapid  non-uniform 
heating  across  the  model  grid. 

Finally,  the  EUVAC  flux  model,  which  is  based  on  a  daily  proxy  and  therefore 
unresponsive  to  the  rapidly  changing  flux  from  a  flare,  was  discarded  in  favor  of  a 
spectrum  derived  from  the  high-resolution  FISM  model.  FISM  was  re-binned  and 
combined  with  the  GOES  0.04  —  0.5  nm  flux  to  generate  a  new  flare  irradiance  ac¬ 
cording  to  the  Solomon  and  Qian  method.  The  model  was  then  tested,  modeling  the 
ionospheric  response  to  the  X17.2  solar  flare  of  28  October,  2003. 

The  model  was  run  with  and  without  the  high-cadence  solar  flux  from  the 
FISM  model  and  GOES-12  solar  X-ray  flux,  also  sampled  at  one  minute.  Electron 
temperature  and  density  variations  across  the  model  domain  were  compared  to  the 
model  run  without  flare  irradiance.  Increases  in  total  electron  content  (TEC),  electron 
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density  at  various  altitudes,  including  at  the  F2  peak  ( NmF2 )  were  observed,  as  were 
positive  and  negative  changes  in  the  electron  temperature  Tc.  The  simulation  results 
were  also  compared  to  several  studies  and  the  new  model  was  found  to  reproduce  the 
observed  increase  in  total  electron  content  within  10%,  especially  near  the  subsolar 
point.  The  change  in  electron  temperatures  in  the  mid-latitudes  are  also  within  an 
order  of  magnitude  of  a  1-dimensional  model.  Although  much  work  remains,  the 
new  version  of  SAMI2  will  be  a  useful  addition  to  the  arsenal  of  models  capable  of 
predicting  the  ionospheric  response  to  solar  flares. 

5.2  Suggestions  for  Future  Work 

The  suggestions  for  future  work  can  be  summarized  into  two  main  categories: 
additional  testing  on  flares  and  other  events,  and  addressing  anomalies  and  issues 
with  the  model  which  cropped  up  in  testing.  While  this  effort  was  successful,  the 
new  version  remains  relatively  untested.  The  original  plan  for  this  project  included 
a  second  flare  event,  also  an  X-class  flare,  however  time  did  not  permit  adequate 
analysis  of  the  second  event.  An  updated  non- flare  version  of  SAMI2  was  also  used  in 
a  study  of  ionospheric  conditions  at  Jicamarca,  Peru,  using  incoherent  scatter  radar 
[Hysell  et  al,  2007].  This  version  of  SAMI2  was  designed  for  parallel  processing  based 
on  the  message  passing  interface  (MPI)  standard.  Although  SAMI2’s  results  were  in 
very  good  agreement  overall,  a  constant  correction  was  applied  to  the  photoelectron 
heating  variable  ( Qphe )  hi  order  to  obtain  the  correct  magnitude  [ Huba ,  2007].  The 
version  of  SAMI2  from  this  paper  generated  lower  temperatures,  however  the  model 
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also  generated  anomalous  patches  of  high  temperature  which  were  deemed  not  to 
be  physically  realistic.  Huba  was  able  to  reproduce  the  instabilities  using  the  open 
source  SAMI2;  however,  the  newer,  parallelized  version  of  the  code  used  in  the  Ji- 
camarca  study  did  not  exhibit  the  instability  [Huba, private  communication].  In  the 
end,  a  comparison  to  [ Hysell  et  al. ,  2007]  was  abandoned,  with  the  hope  of  eventually 
applying  the  changes  from  this  paper  to  the  parallel  version  of  SAMI2.  This  leads  to 
my  first  suggestion,  which  is  to  apply  the  improvements  from  this  study  to  the  MPI 
version  of  SAMI2.  The  open  source  SAMI2  is  inferior  to  the  latest  version,  which  also 
incorporates  a  new  version  of  the  Horizontal  Wind  Model  (HWM07)  [ Hysell  et  al., 
2007].  The  fact  that  the  newest  version  appears  more  stable  is  the  main  reason  for 
this.  Also  a  parallelized  version  would  be  ideal  to  take  advantage  of  today’s  multi-core 
computers.  Another  troubling  issue  with  the  model  is  the  tendency  to  overestimate 
TEC  at  high  solar  activity.  Addressing  this  issue  is  also  necessary  before  the  model 
can  be  fully  utilized.  This  issue  was  not  uncovered  early  enough  in  the  project  to 
properly  address,  especially  since  the  delta  between  the  flare  and  no-flare  runs  was 
in  good  agreement  with  published  results  and  the  other  data  appeared  valid.  The 
problem  may  even  be  in  the  external  routine  used  to  calculate  TEC  from  the  model 
output. 

In  conclusion,  the  new,  improved  version  of  SAMI2  provides  a  new  tool  for 
studying  the  ionospheric  response  to  solar  flares.  The  low  computational  cost  and 
easy  modification  of  the  code  to  add  ion  species,  more  chemistry,  different  wind  or 
electric  fields,  etc.  make  it  an  attractive  research  tool.  Now  it  can  model  the  response 
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to  solar  flares.  It  is  hoped  that  this  additional  capability  will  be  used  within  the 
Department  of  Defense  and  United  States  Air  Force  to  explore  and  better  understand 
the  complex  reaction  of  the  ionosphere  to  a  flare  event.  Understanding  the  inter- 
hemispheric  response  to  flares  is  a  crucial  step  in  predicting  and  mitigating  the  impact 
to  our  communications,  navigations  and  space-borne  assets. 


Appendix  A.  SAMI2  Input  Parameters 


This  explains  and  lists  the  inputs  to  the  modified  version  of  SAMI2  that  were 
used  in  modeling  the  response  to  the  28  October,  2003  solar  flare.  Additional 
parameters  are  listed  and  explained  in  the  README  file  that  accompanies  the  code 
package. 


Table  A.l.  SAMI2  Input  parameters 
for  28  October  2003  solar  flare  simula¬ 
tion. 


Variable 

Value 

USEPE 

.true. 

Turns  photoclectrons  on  and  off 

FLARE 

.true. 

Turns  flare  routine  on  and  off 

fejer 

.true. 

Uses  Fejer  and  Scherliess  electric 
field  model  or  a  sinusoidal  model 

rmin 

90.0 

Height  of  lowest  the  field  line 

rrnax 

7000.0 

Height  of  highest  the  field  line 

nf 

228 

Number  of  field  lines 

nz 

301 

Number  of  grid  point 

fl0p7 

270.9 

Daily  FU)7  index 

fbar 

138.7 

81-day  average  of  daily  Fw,7 

ap 

21 

Ap  index  value 

nion2 

7 

Number  of  ions  (1  —  7) 

tvnO 

1. 

Factor  for  neutral  wind 

tvexbO 

1. 

E  x  B  drift  velocity 

gams 

3 

Grid  spacing  along  the  field  line 

gamp 

3 

Spacing  of  the  field  lines 

cqe 

7.e-14 

Constant  for  non-local 
photoelectron  heating 

USEPE  is  a  logical  parameter  that  turns  on  secondary  ionization  by  photoelec¬ 
trons.  FLARE  enables  the  solar  flare  subroutine  which  then  reads  the  high-cadence 
irradiance  data  starting  at  00  UT  on  the  second  day  of  the  model  run.  The  high 
cadence  data  output  is  started  at  a  time  given  by  beginflr  and  normal  data  output 
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cadence  is  resumed  by  endflr.  Both  beginflr  and  endflr  are  given  in  UT  seconds  + 
86400,  in  order  to  account  for  running  the  model  with  a  24-hour  burn-in  period.  This 
also  requires  the  day  variable  to  be  set  for  the  day  prior  to  the  flare  day.  fejer  is 
also  a  logic  variable  that  specifies  use  of  the  Fejer/ Scherliess  model  of  the  electric 
held  [ Scherliess  and  Fejer,  1999]  used  to  calculate  E  x  B  drift.  Setting  this  to  .false. 
forces  SAMI2  to  use  a  simpler  sinusoidal  drift  model  [Huba  et  al,  2002],  tvnO  is  a 
multiplicative  factor  (0—1)  for  the  neutral  winds.  Setting  this  variable  equal  to  0  turns 
off  the  neutral  winds  and  1  uses  the  unaltered  (100%)  winds  from  the  HWM93  model. 
tvexbO  does  the  same  for  the  E  x  B  drift  velocity,  gams  sets  the  grid  spacing  along  the 
magnetic  held  lines.  The  default  setting  is  3  and  uses  increases  the  spacing  between 
points  as  altitude  increases.  A  lower  number  makes  the  spacing  more  uniform,  gamp 
does  the  same,  but  for  the  spacing  between  held  lines  at  the  magnetic  equator.  As 
altitude  increases  the  spacing  increases  and  a  lower  value  makes  the  spacing  more 
uniform,  cqe  is  the  correction  factor  C  used  in  the  photoclectron  heating  equation, 
equation  3.12.  The  normal  value  ranges  from  3-8  x  10-14. 
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